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Abstract 

The  most  precise  relative  positioning  obtained  using  differential  GPS  depends  on 
accurately  determining  carrier-phase  integer  ambiguities.  In  order  to  achieve  high 
precision,  many  current  static  and  kinematic  surveying  algorithms  use  a  floating-point 
solution  until  enough  information  becomes  available  to  fix  the  carrier-phase  ambiguities 
accurately.  However,  in  dynamic  environments  where  many  brief  measurement  outages 
or  cycle  slips  are  possible,  these  algorithms  may  never  gain  enough  information  to  fix  the 
ambiguities  with  the  confidence  required  for  a  fixed-integer  solution. 

A  mew  method  is  presented  that  uses  a  multiple  model  Kalman  filter  to  resolve  the 
carrier-phase  integer  ambiguities.  This  method  starts  with  the  floating-point  results,  yet 
smoothly  and  rapidly  attains  the  precision  of  the  correct  fixed-integer  solution, 
eliminating  the  need  to  decide  when  to  switch  from  the  floating  to  the  fixed-integer 
solution.  This  method  is  based  on  a  theoretically  correct  blending  of  solutions  from 
multiple  filters,  each  of  which  hypothesizes  a  different  ambiguity  set.  This  new 
technique  is  computationally  efficient,  providing  a  robust  navigation  solution  useful  in 
demanding  applications  such  as  precision  landing  and  autonomous  navigation. 

The  new  method  was  evaluated  during  static  ground  and  flight  tests.  Initial  results 
indicate  that  the  new  method  is  capable  of  quickly  resolving  the  carrier-phase 
ambiguities,  transitioning  smoothly  to  the  true-fixed  integer  solution,  and  providing  a 
highly  accurate  (centimeter-level)  navigation  solution.  However,  the  performance  of  the 
new  method  depends  on  the  correct  carrier-phase  ambiguity  set  being  hypothesized  by 
one  of  the  multiple  filters.  Recommendations  for  further  research  are  included. 
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DEVELOPMENT  AND  TESTING  OF  A 


MULTIPLE  FILTER  APPROACH  FOR  PRECISE  DGPS 
POSITIONING  AND  CARRIER-PHASE  AMBIGUITY  RESOLUTION 


1  Introduction 


1.1  Background 

The  Global  Positioning  System  (GPS),  a  satellite  based  radio  navigation  system,  has 
become  the  accepted,  if  not  the  standard,  navigation  system  for  much  of  the  world.  GPS 
is  used  for  many  navigation  applications:  marine,  airborne,  and  land  navigation,  as  well 
as  in  many  surveying  applications.  This  increasing  dependence  on  GPS  has  also  caused  a 
desire  for  an  increased  accuracy  of  the  GPS  solution. 

GPS-determined  position  and  velocity  can  be  very  accurate,  but  the  accuracy  depends 
on  the  type  of  receiver  and  the  type  of  processing.  Today’s  stand-alone  civilian  receivers 
(using  only  Coarse/ Acquisition  (C/A)  code;  see  Section  1.4.2)  have  an  accuracy  of 
approximately  8-10  meters  root-mean  square  (RMS)  while  military  signal  (P(Y);  see 
Section  1.4.2)  capable  receivers  have  a  slightly  better  accuracy  of  approximately  6-8 
meters  RMS. 
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For  many  applications,  the  accuracy  provided  by  stand-alone  GPS  is  not  sufficient. 
To  attain  better  accuracy,  many  types  of  Differential  GPS  (DGPS)  systems  have  been 
developed.  These  DGPS  techniques  do  not  perform  absolute  positioning,  but  instead 
perform  relative  positioning  between  GPS  receivers.  Differential  GPS  techniques  are 
based  on  the  principle  that  many  of  the  errors  in  the  GPS  signal  are  spatially  correlated 
between  receivers.  Relative  positioning  removes  the  majority  of  the  correlated  errors, 
providing  a  significant  performance  increase  over  stand-alone  GPS.  Differential  GPS 
techniques  can  be  categorized  into  two  main  types  based  on  the  types  of  measurements 
used:  code-based  and  carrier-based.  Representative  position  solutions  and  the 
approximate  accuracy  of  various  types  of  GPS  and  DGPS  systems  are  presented  in  Figure 
1-1  and  Table  1-1.  For  more  detail  the  reader  is  referred  to  references  [20, 34,  35]. 


Time  (Minutes) 


Figure  1-1.  Typical  GPS/DGPS  Solution  Time  Histories 
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Table  1-1.  Typical  GPS/DGPS  Accuracies 


Stand-alone 

Code  based 

Floating-Point 

Fixed-Integer 

GPS 

DGPS 

CPDGPS 

CPDGPS 

6-10  m 

1-3  m 

0.2-0.5  m 

0.01-0.03  m 

Code-based  DGPS  techniques  use  the  GPS  “pseudorange”  measurement.  Code-based 
DGPS  removes  certain  errors  from  the  pseudorange  measurements  to  provide  an  accurate 
position  solution.  The  pseudorange  measurement  is  obtained  by  the  GPS  receiver 
locking  onto  a  GPS  satellite’s  pseudorandom  code.  The  GPS  pseudorange  measurement 
is  essentially  the  difference  between  the  time  of  transmission  of  the  GPS  signal  from  the 
satellite  and  the  time  of  reception  by  the  GPS  receiver.  The  pseudorange  measurement  is 
an  absolute  range  measurement  between  the  receiver  and  the  satellite,  and  therefore 
DGPS  systems  are  relatively  easy  to  implement.  Code-based  DGPS  systems  are  widely  in 
use  today  and  have  a  relative  positioning  accuracy  of  approximately  of  1-3  meters  RMS, 
as  shown  in  Figure  1-1.  The  United  States  Coast  Guard  operates  a  code-based  DGPS 
system  for  precision  marine  navigation  on  most  of  the  country’s  major  waterways  and 
harbors.  Companies  are  also  developing  code-based  DGPS  system  to  monitor  vehicles 
for  a  variety  of  applications  such  as  mining,  farming  and  cross-country  trucking. 
Although  widely  in  use  today,  code-based  DGPS  techniques  are  not  the  focus  of  this 
thesis. 

Although  code-based  DGPS  system  provide  a  significant  increase  in  accuracy  over 
stand-alone  GPS,  carrier-phase  based  DGPS  techniques  must  be  employed  to  achieve  the 
highest  level  of  accuracy  (on  the  order  of  centimeters)  demanded  by  some  applications. 
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Carrier-phase  based  DGPS  techniques  are  based  on  the  phase,  or  carrier  cycles,  of  the 
GPS  carrier  frequency.  A  GPS  carrier-phase  measurement  is  the  difference  between  the 
receiver’s  internally  generated  carrier  signal,  and  the  carrier  signal  received  from  a 
satellite.  Unfortunately,  a  GPS  receiver  cannot  tell  the  difference  between  one  carrier 
cycle  and  the  next.  The  receiver  only  counts  the  change  in  cycles  between  itself  and  the 
satellite.  This  means  the  initial  number  of  cycles,  which  is  by  nature  an  integer  value, 
between  the  receiver  and  satellite  is  unknown,  or  ambiguous.  Because  the  initial  number 
of  cycles  is  unknown,  carrier-phase  measurements  are  not  absolute  range  measurements. 
In  order  to  use  carrier-phase  measurements  for  positioning,  the  carrier-phase  integer 
ambiguities  must  be  somehow  estimated  or  resolved. 

Floating-point  carrier-phase  DGPS  techniques  do  not  resolve  carrier-phase 
ambiguities  to  their  integer  values,  but  simply  estimate  the  carrier-phase  ambiguities, 
generally  using  pseudorange  measurements.  The  carrier-phase  floating-point  ambiguity 
estimates  (which  are  not  integer  values)  are  then  used  to  correct  the  carrier-phase 
measurements  into  absolute  range  measurements.  Floating-point  carrier-phase  DGPS 
techniques  use  both  pseudorange  and  carrier-phase  measurement  to  generate  an  accurate 
relative  position  solution.  Floating-point  carrier-phase  DGPS  systems  have  a  positioning 
accuracy  of  approximately  0.2  -  0.5  meters  RMS.  In  addition,  carrier-phase  solutions 
tend  to  be  “smoother”  than  code-based  solutions,  as  shown  in  Figure  1-1. 

To  obtain  the  highest  level  of  accuracy  available  with  GPS  measurements,  the  carrier- 
phase  integer  ambiguities  must  be  correctly  resolved.  Once  resolved  these  “fixed- 
integer”  based  carrier-phase  DGPS  systems  have  a  positioning  accuracy  of  approximately 
0.01  -  0.03  meters  RMS,  as  shown  in  Figure  1-1.  Carrier-phase  DGPS  systems  are 
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widely  in  use  today  with  surveying  companies.  In  addition,  these  types  of  systems  are 
being  developed  for  precision  navigation,  control,  and  autonomous  landing  of 
commercial  and  military  aircraft. 

Carrier-phase  ambiguity  resolution  generally  consists  of  two  main  functions, 
identifying  possible  ambiguity  sets  and  determining  which  among  these  sets  is  correct. 
Numerous  efficient  methods  have  been  developed  to  address  the  first  part  of  the 
ambiguity  resolution  problem,  and  it  is  not  the  focus  of  this  research.  Additionally,  many 
methods  are  currently  available  to  determine  the  correct  ambiguity  set.  Most  current 
methods  simply  wait  until  enough  information  is  obtained,  by  statistically  comparing 
possible  solutions,  until  the  correct  ambiguity  set  becomes  evident.  The  disadvantage  of 
this  approach  is  it  ignores  the  possibility  of  attaining  a  better  position  solution  until  the 
ambiguities  are  declared  fixed.  Some  new  efforts  have  explored  the  possibility  of  simply 
using  the  best  set,  namely  the  rounded  floating-point  ambiguity  estimate,  at  each  time 
step  to  obtain  the  optimal  solution  [8].  The  disadvantage  of  this  approach  is  that  it  relies 
only  on  a  floating  filter  and  ignores  information  available  from  other  ambiguity  set 
solutions  and  the  information  available  from  the  time  history  of  candidate  ambiguity  sets. 

The  disadvantages  of  both  these  approaches  lend  themselves  to  the  use  of  Kalman 
filter  techniques  in  determining  the  correct  ambiguity  set.  This  work  provides  a 
development  and  analysis  of  using  multiple  model  adaptive  Kalman  filter  techniques 
[22,  25]  for  the  ambiguity  resolution  problem. 
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1.2  Problem  Definition 

The  primary  goal  of  this  thesis  is  the  development  and  testing  of  a  new  methodology 
to  resolve  carrier-phase  ambiguities  using  advanced  Kalman  filtering  techniques.  The 
specific  techniques  are  discussed  in  detail  in  Section  2.5.3.  The  secondary  goal  is  the 
development  and  testing  of  a  carrier-phase  based  DGPS  system  accurate  enough  to  be 
used  as  a  navigation  truth  source  for  airborne  applications.  The  new  methodology  must: 

•  Accurately  determine  the  carrier-phase  ambiguities. 

•  Determine  the  ambiguities  quickly,  or  with  limited  data. 

•  Provide  a  highly  accurate  navigation  solution  at  all  times. 

1.3  Scope 

This  thesis  covers  the  development  and  testing  (both  ground  test  and  flight  test)  of  a 
new  algorithm  designed  to  resolve  carrier-phase  ambiguities.  Performance  evaluations  of 
the  new  algorithm  are  based  on  its  ability  to  resolve  the  carrier-phase  ambiguities  and  its 
position  and  velocity  accuracy  as  compared  to  a  nominal  “truth”  source. 

Software  development  was  conducted  using  MATLAB®  5.3  software  [23]. 
Simulated  GPS  measurements  were  used  to  develop  and  debug  the  software.  Analysis  of 
the  new  methodology  was  conducted  using  actual  GPS  signals.  Ground  testing  consisted 
of  numerous  static  data  sets  collected  using  Ashtech  GPS  Z-Surveyor  receivers  for  both 
the  remote  and  reference  receivers.  Ground  test  evaluations  are  based  on  Ashtech  Office 
Suite  Software-generated  truth  data  [1].  Truth  data  includes  positioning  information  and 
the  correct  carrier-phase  ambiguities.  Flight  test  data  consists  of  numerous  data  sets 
collected  by  flying  a  T-38  aircraft,  again  using  Ashtech  GPS  Z-Surveyor  receivers  [2]  for 
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both  the  reference  and  remote  (airborne)  receivers.  Flight  test  evaluations  are  again 
based  on  Ashtech  Office  Suite  SofMare  generated  truth  data.  Flight  test  truth  data 
includes  positioning  information  and  the  correct  carrier-phase  ambiguities. 

Detailed  information  on  the  data  sets  used  in  this  research  is  presented  in  Chapter  4. 
A  description  of  the  flight  test  profiles  is  presented  in  Appendix  A.  For  more  information 
on  Ashtech  Office  Suite  software,  the  reader  is  referred  to  reference  [1]. 

1.4  Related  Research 

Much  of  the  work  in  this  thesis  is  based  on  current  research  found  in  industry  and 
academia.  This  research  can  be  broken  down  into  two  main  areas:  Kalman  filtering  and 
GPS  (including  various  types  of  Differential  GPS). 

1.4.1  Kalman  Filtering 

Developed  in  1960  by  Kalman,  a  linear  Kalman  filter  is  an  optimal  recursive  data 
processing  algorithm  [24,  25].  All  Kalman  filters  are  based  on  a  model  of  the  system  of 
interest.  The  Kalman  filter  uses  information  about  the  model  of  the  system  as  well  as  the 
measurements  of  the  system  to  generate  optimal  state  estimates.  The  filter  is  optimal 
with  respect  to  virtually  any  criterion,  one  aspect  being  that  the  filter  incorporates  all 
available  information.  The  filter  is  recursive  in  that  the  filter  does  not  require  all  previous 
data  to  be  reprocessed  every  time  a  new  measurement  is  received.  Since  the  introduction 
of  the  linear  Kalman  filter,  new  designs  have  been  developed  to  handle  nonlinear  models 
and  parallel  filter  applications  [22,  25]. 
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1. 4.1.1  Multiple  Model  Adaptive  Estimation 

In  1965,  Magill  [22]  presented  what  is  now  called  a  Multiple  Model  Adaptive 
Estimator  (MMAE).  An  MMAE  is  a  bank,  or  set,  of  individual  Kalman  filters  each  based 
on  a  particular  assumption  about  the  system  of  interest.  Each  Kalman  filter  within  the 
MMAE  is  called  an  “elemental  filter.”  Each  elemental  filter  generates  its  own  optimal 
estimate  of  the  variables  of  interest,  based  not  only  on  the  system  measurements,  but  also 
on  its  own  assumptions  (or  model)  of  the  system.  The  characteristics  of  the  residual 
vector  outputs  of  the  elemental  filters  are  used  to  determine  which  elemental  filter’s 
system  model  best  matches  the  true  system  model  [12].  The  first  practical  MMAE  was 
implemented  in  1977  by  Athens  et.  al.  [3]. 

Optimality  conditions  for  the  MMAE,  similar  to  those  of  a  linear  Kalman  filter,  were 
proven  by  Chang  and  Athans  [5],  but  only  if  the  true  system  exactly  matched  the 
assumption(s)  modeled  by  one  of  the  elemental  filters.  If  the  modeled  system  does  not 
exactly  match  the  real  system,  research  indicates  the  MMAE  will  converge  to  the 
elemental  filter  based  on  “closest  to  true”  system  model,  in  the  Barum’s  distance-measure 
sense  [4,  25]. 

The  Air  Force  Institute  of  Technology  (AFIT)  has  been  a  major  contributor  in  the 
area  of  MMAE  research.  AFIT  has  researched  and  applied  MMAE  algorithms  in  various 
areas  including:  detection  of  failed  sensors  [10,  30],  adaptive  target  tracking  [29], 
adaptive  control  [39],  and  GPS  spoofing  and  jamming  [43]. 

In  1985,  Maybeck  and  Suizu  [28]  addressed  the  problem  of  “beta  dominance”  [3]  in  a 
MMAE.  Beta  dominance  is  the  tendency  of  the  MMAE  to  weight  an  elemental  filter 
incorrectly  due  to  erroneous  residual  contributions.  The  method  developed  by  Maybeck 
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and  Suizu  sets  the  necessary  scalar  coefficient  to  unity;  see  Equation  (2-56).  This  method 
greatly  improves  the  performance  of  the  MMAE  in  certain  applications.  Additional  work 
on  the  beta  dominance  problem  known  as  was  done  by  Menke  [30]. 

In  1987,  Maybeck  and  Hentz  [27]  reduced  the  number  of  necessary  elemental  filters 
by  using  the  concept  of  a  moving  bank  of  filters.  Instead  of  using  a  large  number  of 
elemental  filters,  a  smaller  number  of  points  upon  which  to  base  elemental  filters  are 
moved  as  necessary  through  the  appropriate  parameter  space.  Additional  moving  bank 
MMAE  work  was  done  in  1998  by  Vasquez  [44].  Continuing  the  idea  was  the 
development  of  “  filter  spawning”  in  1998  by  Fisher  [12].  Filter  spawning  involves  using 
a  minimal  set  of  filters  until  it  becomes  necessary  to  use  more  filters  to  define  the 
parameter  space  coverage  better  [12]. 

Recent  work  with  MMAE  design  directly  related  to  GPS  started  in  1992  with 
Vasquez  [43].  The  preliminary  work  used  an  MMAE  structure  for  detection  of  spoofing, 
jamming  or  failures  in  the  GPS  system.  In  1998,  Vanek  [42]  used  an  MMAE  structure 
for  GPS  signal  offset  detection  and  noise  strength  estimation. 

As  shown  above,  there  has  been  an  abundance  of  research  in  the  area  of  MMAE 
design  and  application,  but  an  MMAE  structure  has  never  been  used  to  resolve  carrier- 
phase  ambiguities.  This  thesis  uses  an  MMAE  to  resolve  the  carrier-phase  ambiguities. 

1.4.2  Global  Positioning  System 

The  Global  Positioning  System  (GPS)  is  a  constellation  of  27  radio  navigation 
satellites  (24  primary  plus  three  spares).  Each  satellite  continuously  transmits  a  stream  of 
navigation  codes  and  ranging  information  providing  users  with  data  capable  of  producing 
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the  position,  velocity,  and  time  at  the  user’s  location.  Currently,  most  GPS  receivers 
process  data  from  all  visible  satellites  (up  to  the  number  of  channels  in  the  receiver).  For 
more  details  on  GPS,  the  reader  is  referred  to  references  [20, 34, 35]. 

1.4.2. 1  Carrier-Phase  Ambiguity  Resolution 

Resolving  the  carrier-phase  integer  ambiguities  has  been  an  area  of  research  for  many 
years.  The  simplest  technique  is  to  use  code  measurements  to  determine  the  carrier-phase 
ambiguities  directly.  However,  due  to  the  large  errors  associated  with  code 
measurements,  this  technique  is  usually  not  adequate.  Another  early  attempt  to  resolve 
the  carrier-phase  ambiguities  was  the  Ambiguity  Function  Method  [7].  The  Ambiguity 
Function  Method  is  a  pure  position  domain  search,  and  is  computationally  intensive  and 
not  commonly  used. 

The  Least  Squares  Ambiguity  Search  Technique  (LSAST)  was  developed  by  Hatch 
[17]  in  1990.  In  the  LSAST  method,  a  given  number  of  ambiguity  sets  are  formed  based 
on  the  size  of  the  ambiguity  search  space.  Then  each  ^biguity  set  is  tested  by 
calculating  a  least-squares  position  estimate  and  residual.  The  integer  ambiguity  set  with 
residuals  (the  difference  between  the  true  measurement  and  the  predicted  measurement) 
significantly  lower  than  the  other  ambiguity  sets  is  assumed  to  be  the  correct  set.  This 
method  is  also  computationally  intensive,  but  is  more  efficient  than  the  Ambiguity 
Function  Method.  This  method  is  also  not  commonly  used. 

Other  techniques  use  information  about  the  ambiguities  contained  in  a  covariance 
matrix  to  help  determine  the  correct  ambiguity  set.  The  Fast  Ambiguity  Search  Filter 
(FASF)  was  developed  by  Chen  and  Lapachelle  in  1993  [6].  The  FASF  method  starts 
with  the  floating-point  ambiguity  estimates  and  covariance.  The  method  next  assumes 
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the  best  known,  smallest  covariance,  ambiguity  is  correct  and  rounds  it  to  the  closest 
integer.  Then  new  conditional  covariance  and  ambiguity  estimates  are  developed  for  the 
other  ambiguities  based  on  the  first  ambiguity  being  “correct.”  The  process  is  repeated 
until  a  valid  ambiguity  set  is  found  or  an  integer  ambiguity  is  not  found  within  the  search 
space.  This  results  in  an  efficient  recursive  carrier-phase  ambiguity  resolution  algorithm. 

One  of  the  most  efficient  techniques  is  known  as  the  Least-squares  AMBiguity 
Decorrelation  Adjustment  method  (LAMBDA)  [40,  9].  The  LAMBDA  method  consists 
of  two  steps.  First,  the  ambiguities  are  decorrelated  using  a  linear,  volume-preserving 
integer  transform  known  as  a  “Z-transform”.  Then  a  sequential  conditional  least-squares 
search  is  performed  on  the  transformed  ambiguities.  The  LAMBDA  method  was 
improved  by  Li  and  Gao  [21],  with  the  ability  to  construct  high  dimensional  “Z- 
transforms.” 

Some  recent  research  has  been  performed  in  applying  advanced  Kalman  filtering 
techniques  to  the  ambiguity  resolution  problem.  Cox  and  Brading  [8]  proposed  using  the 
outputs  of  a  Kalman  filter  to  determine  the  optimal  integer  ambiguity  set  and  the  optimal 
navigation  solution.  This  approach  considers  the  integer  values  closest  to  the  Kalman 
filter  estimates  as  the  “the  optimal”  ambiguity  set  and  does  not  use  MMAE  techniques  to 
determine  the  correct  ambiguity  set.  Finally,  in  2000  Williamson  and  Speyer  [45]  used 
an  MMAE  structure  to  detect  carrier-phase  cycle  slips  in  an  integrated  GPS/INS  (Inertial 
Navigation  System),  but  not  to  resolve  the  integer  ambiguities. 
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1.5  Methodology 

The  research  presented  in  this  thesis  begins  with  Kalman  filtering  and  MMAE 
designs  and  applications  presented  in  Maybeck  [27,  28],  Magill  [22],  Vasquez  [44], 
Miller  [31],  Vanek  [42]  and  Fisher  [12].  This  produced  the  basic  MMAE  design  used  in 
this  thesis.  In  addition,  a  detailed  study  was  conducted  on  GPS,  DGPS  and  carrier-phase 
ambiguity  resolution  techniques  presented  by  Kaplan  [20],  Raquet  [36,  37],  Hatch  [17], 
Chen  and  Lapachelle  [6],  and  Teunissen  [40]. 

A  basic  DGPS  Kalman  filter  was  implemented  in  MATLAB®  [23]  using  a  design 
similar  to  the  one  outlined  by  Raquet  [36].  The  filter  was  expanded  to  use  single  or  dual 
frequency  data,  and  other  measurement  combinations  as  specified  by  the  user. 

The  development  of  the  overall  algorithm  began  after  verifying  the  performance  of 
the  basic  DGPS  model.  The  MMAE  model  was  designed  using  the  basic  MMAE 
structure  found  in  Maybeck  [25].  It  contains  up  to  100  elemental  filters,  each  based  on  a 
different  ambiguity  set.  Ambiguity  set  generation  for  this  thesis  was  developed  using  a 
combination  of  the  LAMBDA  method  “Z-transform”  and  an  FASF  subroutine.  The 
combination  of  these  two  methods  is  detailed  in  Section  2. 5.3.1. 

All  analysis  of  the  new  methodology  was  conducted  using  actual  GPS  signals. 
Ground  tests  were  conducted  with  Ashtech  GPS  Z-Surveyor  receivers  for  both  the  remote 
and  reference  receivers.  The  performance  of  the  algorithm  was  analyzed  based  on  four 
different  ground  test  cases,  with  baselines  varying  from  250  meters  to  100  kilometers. 
Ground  test  evaluations  are  based  on  Ashtech  Office  Suite  Software-generated  “truth” 
data,  which  includes  positioning  information  and  the  correct  carrier-phase  ambiguities. 
Flight  tests  were  conducted  by  flying  a  T-38A  aircraft,  also  with  Ashtech  GPS  Z- 
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Surveyor  receivers  for  both  the  remote  and  reference  receivers.  The  performance  of  the 
algorithm  was  analyzed  based  on  numerous  different  flight  test  cases,  with  baselines 
varying  from  10  nautical  miles  to  100  nautical  miles  at  altitudes  from  10,000  to  30,000 
feet.  Flight  test  evaluations  are  based  on  GPS  Aided  INS  Navigation  Reference 
(GAINR)  generated  truth  data,  which  includes  position  and  velocity  information,  and 
Ashtech  Office  Suite  Software  truth  data  containing  the  correct  carrier-phase  ambiguities. 
Finally,  analysis  and  plotting  routines  were  built  in  MATLAB®  for  a  complete  tabular 
and  graphical  presentation  that  enabled  comparisons  of  the  algorithm  to  known  truth 
sources,  and  for  side-by-side  case  comparisons. 

1.6  Thesis  Overview 

Chapter  2  presents  the  required  Kalman  Filter,  MMAE,  and  GPS  theory  used  to 
develop  the  mathematical  models  for  this  research.  Kalman  filter  and  MMAE  theory  are 
defined  and  presented  as  is  the  theory,  mathematical  relationships,  and  measurement 
types  for  GPS.  Finally,  the  concepts  and  methods  used  for  candidate  ambiguity  set 
generation,  the  LAMBDA  and  FASF  methods  are  described. 

Chapter  3  describes  the  flow  of  the  new  algorithm.  Detailed  model  and  measurement 
equations  are  presented  for  the  DGPS  filter  and  the  MMAE  design.  The  method  used  to 
generate  candidate  ambiguity  sets  is  described  in  detail,  as  are  the  calculations  required  to 
determine  the  correct  ambiguity  set.  Finally,  concepts  necessary  to  understand  the  output 
of  the  overall  algorithm  are  presented. 
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Chapter  4  presents  a  performance  analysis  of  each  test  case  in  this  research.  Position 
accuracy  analysis  and  ambiguity  resolution  analysis  is  presented  for  ground  and  flight 
tests. 

Chapter  5  presents  a  summary  of  this  research  effort.  In  addition,  conclusions  and 
recommendations  for  future  research  are  presented. 
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2  Theory 


2.1  Overview 

This  chapter  describes  the  linear  Kalman  filter,  the  Extended  Kalman  Filter  (EKF) 
and  the  Multiple  Model  Adaptive  Estimation  (MMAE)  structure  in  more  detail.  In 
addition,  this  chapter  develops  theory  behind  GPS,  Differential  GPS  (DGPS)  and  Carrier- 
Phase  Differential  GPS  (CPDGPS). 

2.2  Kalman  Filters 

A  linear  Kalman  filter  is  an  optimal  recursive  data  processing  algorithm  [24,  25]. 
The  filter  is  optimal  with  respect  to  virtually  any  criterion,  one  aspect  being  the  filter 
incorporates  all  available  information.  The  filter  is  recursive  in  that  the  filter  computes  a 
new  estimate  at  discrete  time  intervals  whenever  a  new  measurement  is  available.  It  does 
not  need  all  previous  measurements  to  produce  an  optimal  state  estimate,  as  would  a 
batch  processing  algorithm.  In  addition,  the  Kalman  filter  is  a  data  processing  algorithm, 
meaning  it  is  simply  a  computer  program  processing  discrete-time  measurement  samples. 

Every  Kalman  filter  is  based  on  a  model  of  the  system  of  interest.  The  Kalman  filter 
uses  system  measurements  and  the  system  model  to  generate  optimal  state  estimates.  The 
Kalman  filter  uses  all  available  measurements,  regardless  of  their  accuracy,  to  generate 
an  optimal  estimate  of  the  variables  of  interest  based  on  knowledge  of  the  system 
dynamics  and  measurement  characteristics,  the  statistical  description  of  the  system 
noises,  measurement  errors,  and  model  uncertainties  [46]. 
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The  Kalman  filter  can  be  described  by  two  stages  of  operation:  the  propagation  stage, 
and  the  update  stage.  During  the  propagation  stage,  the  filter  computes  a  prediction  of 
the  system  states  at  some  future  time  based  on  its  own  imperfect  dynamics  model  of  the 
system.  During  the  update  stage,  the  filter  uses  noise-corrupted  measurements  to  refine 
its  predicted  state  estimates  further,  generating  the  optimal  state  estimate.  A  complete 
derivation  can  be  found  in  Stochastic  Models,  Estimation,  and  Control,  Volume  1  [24]. 

2.2.1  State  and  Measurement  Model  Equations 

Following  the  development  discussed  in  references  [24,  25,  31],  assume  that  the 
system  dynamics  can  be  modeled  as  a  linear,  time-invariant  system.  The  state  equation 
takes  the  form: 

x(0  =  ¥{t)x(t)  +  B(0u(0  +  G(r)w(0  (2-1) 

where 

x(0  =  the  n-dimensional  system  state  vector 

F(0  =  the  n-hy-n  state  dynamics  matrix 

B(0  =  the  n-by-r  control  input  matrix 

u(r)  =  the  r-dimensional  control  input  vector 

G(r)  =  the  n-hy-s  noise  input  matrix 

w(0  =  the  .s-dimensional  dynamics  driving  noise  vector 

The  dynamics  driving  noise  vector  w(0  is  assumed  to  be  white  Gaussian  noise  with 
mean  and  covariance  defined  by: 
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where  t  has  units  of  time. 


£{w(0}=  0 

£{w(r)w^  (t  +  t)}=  Q(t)S{T) 


(2-2) 

(2-3) 


Next,  assume  the  measurements  can  be  modeled  as  a  linear,  time-invariant,  discrete 
time  system.  The  measurement  equation  takes  the  form: 

z(r,. )  =  H(r,.  )x(r,. )  +  v(r,. )  (2-4) 

where 


z(r,. )  =  the  /n-dimensional  measurement  vector 

H(r, )  =  the  m-by-n  system  output  matrix 

v(r, )  =  the  m-dimensional  discrete  time  noise  input  vector 

The  discrete  time  noise  input  vector  v(r.)  is  also  assumed  to  be  white  Gaussian 
noise.  The  mean  and  covariance  are  defined  by: 

£{v(r,)}=0  (2-5) 


R(r,. )  for  t.  =  tj 
0  for  r.  ^  tj 


(2-6) 


Noises  w(0  and  v(r, )  are  assumed  independent.  Thus: 


4w(0v^(f,)}=0 


(2-7) 


for  all  t  and  r,. 


2.2.2  Kalman  Filter  Equations 

The  Kalman  filter  is  propagated  forward  from  time  U.i  to  t,  starting  from  initial 
conditions  x(ro )  and  ?(?(, ) ,  or  from  the  last  update  cycle  state  and  covariance  estimates. 
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The  Kalman  filter  is  propagated  forward  by: 

i(t/r,_,)  =  F(0x(t/t,_,)  +  B(0u(0  (2-8) 

P(t/t,_,)  =  F(0P(t/t,..i)  +  P(t/?M)F"(0  +  G(t)Q(0G^(0  (2-9) 

Before  being  used  in  the  update  stage  of  the  Kalman  filter  the  state  and  covariance 
differential  equations  are  solved  by: 

x(tr)  =  <D(r,,t,.,)x(t;.i)  +  f  0(t,,T)B(T)u(T>^T  (2-10) 

P(t:)  =  0(f,,,t,_,)P(C,)<D"(t,.,f,_i)  +  f  a>(t,,r)G(T)Q(r)G"(T>D"(t,,T)(!>r  (2-11) 

•"i-i 

where  <&(?,.,?,_,)  is  the  state  transition  matrix  from  time  tui  to  t,-  and  is  calculated  by 
=  <I)(A0  =  where  At  =  t,.  . 

In  Equations  (2-10)  and  (2-11),  the  -  and  +  superscripts  indicate  times  immediately 
before  and  after  measurements  are  incorporated,  respectively. 

The  Kalman  filter  is  updated  by: 


A(t,)  =  H(t,)P(t:)H"(t,)  +  R(t,) 

(2-12) 

K(t,)  =  P(tr)H^(t,)A(t,)-' 

(2-13) 

r(t,)  =  z,-H(t,)x(t:) 

(2-14) 

x(t;)  =  x(t:)  +  K(t,,)r(t,.) 

(2-15) 

P(t;)  =  P(t-)-K(t,)H(t,)P(t:) 

(2-16) 

In  a  properly  designed  filter  the  residual  vector  r(t,.)  is  zero  mean  with  a  filter- 
predicted  covariance  of  A(t,.) .  The  primary  outputs  of  a  Kalman  filter  are  the  optimal 
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state  estimates  x{tl' )  and  their  filter  computed  covariance  P(t,^) .  In  MMAE  designs,  the 
residual  vector  r(i, )  and  predicted  covariance  A(0  are  also  useful. 

2.3  Extended  Kalman  Filters 

A  linear  Kalman  filter  is  only  useful  when  both  the  dynamics  and  measurement 
equations  of  a  system  are  a  linear  combination  of  the  states,  known  control  inputs  and 
noises  modeled.  However,  the  equations  relating  the  state  and  measurement  variables  for 
most  real  system  (including  the  GPS  systems  in  this  thesis)  are  nonlinear.  To  handle 
these  nonlinearities,  the  Extended  Kalman  Filter  (EKF)  is  used. 

The  concept  behind  an  EKF  is  to  relinearize  the  nonlinear  model  about  the  current 
best  estimate  x(t^ ) .  This  relinearization  about  each  new  estimate  allows  the  use  of  linear 

assumptions  and  equations.  The  EKF  is  a  first-order  approximation  to  an  optimal 
nonlinear  filter.  It  is  implemented  in  a  very  similar  manner  as  the  linear  Kalman  filter 
through  a  basic  propagation  stage  and  an  update  stage.  A  complete  derivation  of  the 
EKF  can  be  found  in  Stochastic  Models,  Estimation,  and  Control,  Volume  II  [25]. 


2.3.1  State  and  Measurement  Model  Equations 

Following  the  development  discussed  in  reference  [25],  assume  that  the  system 
dynamics  can  be  modeled  as  a  nonlinear  state  equation,  which  takes  the  form: 

x(t)  =  f[x(t),u(0,t]  +  G(t)w(r)  (2-17) 

where 

x(0  =  the  n-dimensional  state  dynamics  vector 
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f[x(if)^u(f),?]  =  the  state  dynamics  vector,  which  is  a  linear  or  nonlinear  function  of 

the  n-dimensional  state  vector  x(t),  the  r-dimensional  control  input 

function  u(0  and  of  continuous-time  t. 

G(t)  =  the  n-hy-s  noise  distribution  matrix 

w(if)  =  the  i'-dimensional  dynamics  driving  noise  vector 

The  dynamics  driving  noise  vector  is  assumed  to  be  a  white  Gaussian  noise  process 

whose  mean  and  covariance  are  defined  by: 

£:{w(r)}=  0  (2-18) 

(t  +  r)}=  Q(t)S{T)  (2-19) 

where  rhas  units  of  time. 

Next,  assume  the  EKF  discrete  time  measurements  can  be  modeled  as  a  nonlinear 
state  equation.  It  takes  the  form: 

z(t,. )  =  h[x(t. )  t.  ]  +  \{t, )  (2-20) 

where 

z(fi )  =  the  m-dimensional  measurement  vector  at  time  t, 

h[x(t,)t,.]=  a  known  linear  or  nonlinear  m-dimensional  function  of  x(f,. )  and  U 

\(t. )  =  the  m-dimensional  discrete-time  noise  input  vector 


The  discrete-time  noise  input  vector  v(t,.)  is  also  assumed  to  be  white  Gaussian 
noise.  The  mean  and  covariance  are  defined  by: 


£{y(O}=0 

£{v(Ov' (>,)}=  15'^ 

LV,/  KjJi  ^0  for  t,  r. 


(2-21) 

(2-22) 
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2.3.2  State  and  Measurement  Model  Linearization 


Assuming  Equations  (2-17)  and  (2-20)  are  nonlinear,  they  must  be  linearized  in  order 
to  produce  an  optimal  estimate  of  the  state  vector  x(t),  to  first  order.  The  linearization 
technique  used  in  reference  [25]  involves  perturbations  about  a  nominal  state  trajectory  to 
produce  a  first  order  linear  approximation  of  Equations  (2-17)  and  (2-20)  for  use  in  an 
EKE  [46]. 

Assume  a  nominal  state  trajectory  x„  (•)  can  be  generated  from  the  initial  conditions 
(to )  ^no  using  the  following  differential  equation. 

x„(t)  =  f[x„(t)u(0,t]  (2-23) 

In  addition  to  the  nominal  state  trajectory,  assume  a  sequence  of  nominal 
measurements.  They  can  be  defined  by: 

(2-24) 

To  produce  the  perturbation  state  S\(t)  the  nominal  trajectory  (2-23)  is  subtracted 
from  the  system  model.  Equation  (2-17),  to  give: 

Sk{t)  =  [x(t)-x„ (t)]  =f [x(r),u(0,t]  -  f [x„ (t)u(0,t]  +  G(t)w(r)  (2-25) 

A  Taylor  series  expansion  of  f[x(r)u(0,t]  about  x„(0  gives: 


f  [x(t )  u(0,  t]  =  f  [x„  (0,  u(0,  t]  + 


^[x(0,u(0,t] 

Sk 


[x(t)-x„(t)]  +  h.oS. 

X=X„(I) 


(2-26) 


where  “h.o.t.”  represents  the  higher  order  terms  in  the  expansion  with  powers  greater  than 
one.  Finally,  combining  Equations  (2-25)  and  (2-26),  subtracting  like  terms  and  then 
reducing  it  to  a  first  order  approximation  gives: 
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6k{t)  =  ¥[f,  x„  (?)}5x(0  +  G(?)w(0  (2-27) 

Here,  F[';x.W]  is  the  n-hy-n  partial  derivative  matrix  of  f  with  respect  to  x, 


evaluated  along  the  nominal  trajectory  and  defined  by: 


Sk 


(2-28) 


X=X„(0 


This  is  an  acceptable  approximation  to  the  true  solution  provided  the  deviations  from  the 
nominal  solution  are  small  enough  that  the  higher  order  terms  in  the  Taylor  series 
approximation  of  Equation  (2-26)  are  negligible. 

In  a  similar  manner  the  measurement  perturbation  <S!(t,)  of  the  measurement  vector 
z(ti)  is  formed  by  subtracting  Equation  (2-24)  from  (2-20)  to  yield: 

(5z(t, )  =  [z(t,. )-  z„  )]  =  h[x(t, )  t,.  ]  -  h[x„  )  t,.  ]  +  v(t,. )  (2-29) 

A  Taylor  series  expansion  of  h[x(]f,)t,  ]  about  x„(t,  )  gives: 


(2-30) 


Combining  Equations  (2-29)  and  (2-30),  subtracting  like  terms  and  then  reducing  it  to  a 
first  order  approximation  gives: 

<5z(t,)=H[t;;x„(t,.)]^(r,.)+ v(t,.)  (2-31) 

Here,  H[ti;  x„(r,)]  is  the  m-by-n  partial  derivative  matrix  of  h  with  respect  to  x,  evaluated 
along  the  nominal  trajectory  and  defined  by: 


TTr  / 

H[t,.;x„(t,)]= - ^ 


Sk 


(2-32) 
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Again,  this  is  an  acceptable  approximation  to  the  true  solution  provided  the  deviations 
from  the  nominal  solution  are  small  enough  that  the  higher  order  terms  in  the  Taylor 
series  approximation  of  Equation  (2-30)  are  negligible. 

The  system  model  in  Equations  (2-27)  and  (2-31)  is  optimal  to  first  order  due  to  the 
truncated  Taylor  Series,  and  has  been  linearized  to  produce  error  state  equations.  They 
can  now  be  used  to  implement  a  linearized  Kalman  filter  to  perform  all  state  propagation 
and  update  equations.  The  total  state  estimate  can  be  calculated  by  adding  the  nominal 
state  value  to  the  error  state  estimate. 

Recall  that  an  EKF  relinearizes  the  nonlinear  model  about  the  current  best  estimate 
\{t^ )  and  a  trajectory  emanating  from  it,  as  opposed  to  the  nominal  state  just  described. 

This  effectively  incorporates  a  new  and  better  state  trajectory  into  the  filter  model  at  each 
measurement  time,  allowing  the  use  of  linear  perturbation  techniques. 


2.3.3  Extended  Kalman  Filter  Equations 

The  EKF  state  and  covariance  estimates  are  propagated  forward  from  time  t,.]  to  U  by 
integrating  from  the  initial  conditions  x(to )  and  P(t(, ) ,  or  from  the  last  update  cycle 
state  and  covariance  estimates.  The  propagation  equations  are  given  by: 

X(t  /  )  =  f  [x(t  /  ),  u(0,  t]  (2-33) 

P(t/r,_i)  =  F[r;x(t/t,_,)]P(t/t,.,)  +  P(t/t,_,)F^[t;x(t/r,_i)]  +  G(OQ(r)G^^  (2-34) 
where  the  matrix  F[t;  x(t  /  r,._j )]  is  given  by: 
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In  the  above  equations,  the  term  t  / is  the  value  of  a  given  variable  at  time  t. 


conditioned  on  the  measurements  taken  through  timet,._i . 

If  not  at  time  to,  the  initial  conditions  at  the  start  of  the  propagation  cycle  are  defined  by: 

(2-36) 

P{!,.,/(,_,)=Pfc,)  (2-37) 

After  integrating  Equations  (2-33)  and  (2-34)  forward  to  the  next  sample  time,  the 
following  state  and  covariance  values  are  used  in  the  next  measurement  update  cycle: 

x('r)=  x((| /t,.,)  (2-38) 

P('r)=P(',(t,-.)  (2-39) 

The  EKF  is  updated  by: 

K(^)=  p(';)H4,;4r)lHk4r)H- )h4',;4,')]+  R(',)r  (2-40) 

*k)=  4r)+  K((,){z,  (2-41) 

Pk  )  =  p(T  )  -  K(<,  )Ht,  -.xif-  )H- )  (2-42) 


where  the  matrix  H[t,. ;  )J  is  given  by: 

H[.,;xk)]=^_ 


(2-43) 


The  updated  state  estimate  and  covariance  estimates  are  used  as  the  new  initial  conditions 
for  the  propagation  cycle  of  the  sampling  interval: 


x(r,./r,.)=  x(t") 

p(t,/t,)  =  p(t;) 


(2-44) 

(2-45) 
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2.4  Multiple  Model  Adaptive  Estimation 

There  are  many  different  types  of  Multiple  Model  Adaptive  Estimators.  The  MMAE 
used  in  this  research  follows  the  development  of  the  multiple  model  filter  algorithm 
derived  in  reference  [25].  A  block  diagram  of  the  multiple  model  filter  algorithm  is 
shown  below  in  Figure  2-1. 


Figure  2-1.  Multiple  Model  Filtering  Algorithm 


An  MMAE  is  a  bank,  or  set,  of  K  Kalman  filters,  each  designed  using  a  different 
parameter  vector  (denoted  ),  or  model,  of  the  system  of  interest.  Each  one  of  the 

individual  Kalman  filters  is  called  an  elemental  filter.  The  residual  vector,  ,  and  its 
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covariance,  A^ri,.),  are  measures  of  how  well  the  system  modeled  in  the  elemental  filter 
matches  the  true  system.  Therefore,  if  each  elemental  filter  is  given  the  same 
measurements,  then  by  monitoring  the  residuals  it  is  possible  to  determine  which 
elemental  filter’s  model  best  matches  the  true  system. 


2.4.1  Elemental  Filters 

As  stated  above,  each  of  the  K  elemental  filters  in  the  MMAE  are  based  on  a  different 
parameter  vector  a  or  model  of  system  of  interest.  For  example,  the  filter  is  based  on 
the  parameter  vector  .  The  elemental  filters  are  each  propagated  forward  by: 


K ih  ti-i )x, (tU ) + (t,_i )u(t,_i ) 


(2-46) 

(2-47) 


The  elemental  filters  are  each  updated  by: 

A,  (r, )  =  H,  (t,.  )P,  {t:  )H[  (t, )  +  R,  (r, )  (2-48) 

K*(t,)  =  P,(tr)H[(t,)A,(t,)-'  (2-49) 

r*  (ti )  =  z,  -  ih )K  (^r )  (2-50) 

x,(t;)  =  x,(t-)  +  K,(t,)r,(t,)  (2-51) 

P,  {t; )  =  P,  (t : )  -  K ,  (t,.  )H ,  (r,.  )P,  (t : )  (2-52) 

where  the  subscript  k  refers  to  the  elemental  filter. 


It  is  important  to  understand  that  at  every  sample  time  t. ,  each  elemental  filter 
generates  an  optimal  state  estimate  x*.(t,^) ,  covariance  P)fc(t/) .  residual  vector  (t,.) ,  and 
covariance  A^(t,.)  based  on  its  internal  system  model  and  system  measurements. 
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2.4.2  Conditional  Probability  Generation 

The  residuals  of  the  elemental  filters  are  evaluated  by  the  conditional  probability 
computation  shown  in  Figure  2-1.  The  conditional  probability  of  the  elemental  filter 
at  time  t. ,  is  defined  as: 

/?^(t,.)=Prob(a  =  a;t|Z(t,.)  =  Z,.)  (2-53) 

The  conditional  probability  is  the  probability  that  the  random  variable  a ,  representing  the 
system  characteristics,  is  the  realization  of  the  variable  in  the  jt"*  elemental  filter  a^ , 

given  the  time  history  of  measurements  up  to  and  including  that  taken  at  time  r,. .  The 
time  history  of  measurements  is  represented  by: 

Z(t.)  =  [z(tj)^ \z{t^y 

As  shown  in  [25],  the  conditional  probabilities  can  be  computed  as: 

^  f  E((,)|a,Z(r,..,)  (^iK  ’  "^i-OPk  (^i-l ) 

Pk\^i)  r  I 

2^  y=i  f 1(1,  )|a,Z(ii.i )  1^;  ’  ^i-1 )  Pj  (h-l ) 

The  denominator  in  Equation  (2-55)  acts  as  a  scaling  function  ensuring  that  the 
probability  of  each  elemental  filter  is  properly  defined  in  the  sense  that: 

p,  it, )  >  0  for  all  K  and  Pk  ih  )  =  1 

The  conditional  density  function  in  Equation  (2-55)  is  for  the  current  measurement  based 
on  the  elemental  filter’s  parameter  vector  a^t  and  the  time  history  of  measurements.  This 
density  function  is  represented  and  calculated  by: 


(2-54) 


(2-55) 
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1 


T7T*exp{*} 


(2;rr"'|A,(f,)| 

{•}  =  {-^r,^(^,)A,a,r'r,a,) 


(2-56) 


where  m  =  the  number  of  measurements.  Equation  (2-56)  is  the  conditional  density 
function  for  the  residuals  of  the  k  elemental  filter,  conditioned  on  and  the  time  history 
of  measurements  up  to  time  t,-./,  where  r^(t,)  and  are  given  by  Equations  (2-12) 

and  (2-14).  If  designed  correctly,  the  elemental  filter  representing  the  true  system  (i.e., 
having  the  correct  parameter  vector  )  will  have  zero-mean  white  Gaussian  residuals. 


with  values  consistently  more  in  consonance  with  that  filter’s  computed  (t,.) ,  than  the 

residuals  of  other  elemental  filters.  Therefore,  based  on  Equations  (2-55)  and  (2-56),  the 
probability  of  the  correct  filter  will  increase  while  the  probability  of  the  other  filters  will 
decrease. 


2.4.2.1  Lower  Probability  Bounds 

The  conditional  probability  calculation  for  each  elemental  filter  is  recursive  in  nature. 
Notice  in  Equation  (2-55)  that,  if  an  elemental  filter’s  conditional  probability  becomes 
zero,  it  remains  zero  for  all  time.  This  effect  causes  the  MMAE  to  ignore  information 
from  that  particular  elemental  filter.  To  avoid  this  problem,  a  lower  probability  bound  e 
(generally  0.001  <e<  0.01 )  is  set  for  each  elemental  filter.  If  the  calculated  conditional 
probability  for  the  filter  is  less  than  the  lower  bound,  the  filter’s  probability  is  set  to  the 
lower  bound  e  and  then  all  computed  Pi(t,  )’s  are  rescaled  so  that  they  sum  to  one. 

Large  E  values  force  inappropriate  weighting  of  estimates  associated  with  the  wrong 
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filter,  while  very  small  e  values  slow  down  the  response  of  the  MMAE  to  true  parameter 
or  system  changes. 

2.4.3  Overall  MMAE  State  Estimate 

The  output  of  the  multiple  model  filtering  algorithm  is  a  probabilistically  weighted 
average  of  the  outputs  of  the  K  elemental  filters.  The  overall  state  estimate  of  the 
MMAE  is  calculated  by: 

=  (2-57) 

The  conditional  covariance  of  the  MMAE  is  calculated  by: 

''’(^(=SLp*('<)(p*('/')+(**<'r)-*(<*)B.«;)-*(r)ri  (2-58) 

It  is  important  to  understand  that  the  output  of  the  MMAE  is  a  weighted  blending  of  the 
outputs  of  all  of  the  elemental  filters,  including  filters  based  on  the  less  correct 
assumptions. 

2.4.4  MMAE  Convergence 

The  ability  of  the  MMAE  to  converge  to  the  correct  hypothesis  is  guaranteed 
theoretically  under  certain  specified  conditions  [5].  Some  of  these  proofs  require  that  the 
characteristics  a^t  modeled  by  an  elemental  filter  exactly  match  the  true  system 
characteristics.  Additional  research  has  indicated  that,  if  the  modeled  system  does  not 
exactly  match  the  real  system,  the  MMAE  will  converge  to  the  elemental  filter  based  on 
“closest  to  true”  system  model  [25,  26,  27].  However,  none  of  the  work  to  date  gives  any 
guarantee  on  the  convergence  rate  of  the  MMAE. 
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2.5  Global  Positioning  System 

The  Global  Positioning  System  (GPS)  is  a  constellation  of  27  radio  navigation 
satellites  (24  primary  plus  three  spares).  Each  satellite  continuously  transmits  a  stream  of 
navigation  codes  and  ranging  information  providing  users  with  data  capable  of  producing 
the  position,  velocity,  and  time  at  the  user’s  location.  Currently,  most  civilian  GPS 
receivers  process  data  from  all  visible  satellites,  up  to  the  number  of  channels  in  the 
receiver.  The  GPS  development  in  this  work  follows  the  information  provided  in 
references  [11,  20,  36,  37].  For  a  detailed  overview  of  GPS,  the  reader  is  referred  to 
reference  [20] 

2.5.1  The  GPS  Signal 

The  GPS  signal  is  transmitted  continuously  from  all  satellites  on  two  carrier 
frequencies:  an  LI  carrier  (1575.42  MHz)  and  an  L2  carrier  (1227.6  MHz).  Both  the  LI 
and  L2  frequencies  are  modulated  with  specific  pseudorandom  noise  (PRN)  codes  unique 
to  each  satellite.  The  coarse/acquisition  (C/A)  PRN  codes  are  1023  bit  sequences 
modulating  the  carrier  at  1.023  MHz  (repeated  every  millisecond).  The  precision  (P) 
code  is  a  much  longer  code  that  repeats  every  267  days.  Each  GPS  satellite  is  assigned  a 
one-week  segment  of  the  P-code,  which  modulates  the  carrier  at  10.23  MHz.  The  Y-code 
is  an  encrypted  version  of  the  P-code  and  is  referred  to  as  the  P(Y)-code.  The  P(Y)-code 
is  intended  for  military  users  and  allows  for  the  anti-spoof  (AS)  mode  of  operation.  The 
LI  frequency  contains  both  C/A-code  and  a  P(Y)-code,  while  the  L2  frequency  only 
contains  P(Y)-code.  A  50  Hz  navigation  message  is  also  modulated  onto  both  carriers. 
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This  navigation  message  contains  information  for  all  satellites  such  as  satellite  ephemeris, 
atmospheric  propagation,  satellite  clock  bias,  and  almanac  data.  Table  2-1  below 
summarizes  the  GPS  signals  [11]. 


Table  2-1.  GPS  Signal  Summary 


Carrier 

Wavelength 

Modulation 

Frequency 

Chip  Length 

C/A  Code 

1.023  MHz 

293  m 

LI 

1575.42  MHz 

19  cm 

P  Code 

10.23MHz 

29.3  m 

Nav  Message 

50  Hz 

L2 

1227.60  MHz 

24  cm 

P  Code 

10.23MHz 

29.3  m 

Nav  Message 

50  Hz 

GPS  receivers  are  called  single  or  dual  frequency  receivers  based  on  their  ability  to 
receive  either  LI  signal  or  both  LI  and  L2  signals.  Receivers  correlate  the  received  code 
with  an  internally  generated  copy  of  the  appropriate  PRN  code.  Single  frequency 
receivers  use  only  C/A  code  while  dual  frequency  receivers  use  both  C/A  and  P  codes. 
Receivers  that  cannot  decrypt  the  P(Y)  code  signal  use  semi-codeless  techniques  to  get 
range  information  from  the  encr3q)ted  signal  [41]. 

At  the  start  of  this  research,  the  GPS  signal  was  intentionally  degraded  by  a  process 
known  as  selective  availability  (SA).  Selective  availability  was  applied  to  the  signal 
before  being  transmitted  from  the  satellite  and  was  either  a  pseudo-random  dither  of  the 
satellite  GPS  clock,  or  errors  added  into  the  ephemeris  information  broadcast  in  the 
navigation  message.  This  denied  full  position  accuracy  to  users  not  capable  of  decrypting 
the  P(Y)-code.  SA  decreased  positioning  accuracy  to  approximately  100  meters  (2-D 
RMS)  with  SA  from  approximately  10  meters  (2-D  RMS)  without  SA.  In  the  spring  of 
2000,  SA  was  turned  off  and  full  accuracy  is  available  to  all  users;  however  civilian  users 
still  only  have  access  to  C/A  code  signals. 
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2.5.2  GPS  Measurements 


To  calculate  the  user’s  position,  the  GPS  receiver  must  measure  the  range  between 
itself  and  each  satellite  in  view.  The  GPS  signal  contains  two  main  measurement  types: 
code  and  carrier-phase. 

2.5.2. 1  Code  Measurements 

The  GPS  code  measurement  is  called  a  called  a  pseudorange  because  it  is  a 
combination  of  the  true  range  to  the  satellite  and  errors  in  the  GPS  signal.  It  is  essentially 
the  difference  between  the  time  of  transmission  from  the  satellite  of  the  GPS  signal  and 
the  time  of  reception  by  the  GPS  receiver.  The  pseudorange  measurement  is  defined  as: 

yO  =  r  +  +  7’  +  7  m 1)  SA  (2-59) 

where 

p  =  GPS  pseudorange  (m) 
r  =  true  user  to  satellite  range  (m) 
c  =  speed  of  light  (m/s) 

=  user  (receiver)  clock  error  (s) 

=  satellite  clock  error  (s) 

T  =  errors  due  to  tropospheric  delay  (m) 

1  =  errors  due  to  ionospheric  delay  (m) 
m  =  errors  due  to  signal  multipath  (m) 

V  =  errors  due  to  receiver  noise  (m) 

SA  =  errors  sue  to  selective  availability  (m) 
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2.5.2.2  GPS  Carrier-Phase  Measurements 


A  GPS  carrier-phase  measurement  is  the  difference  between  the  receiver’s  internally 
generated  carrier  signal,  and  the  carrier  signal  received  from  a  satellite.  The  carrier-phase 
measurement  (in  cycles)  is  defined  as; 

(/)  =  -(r  +  c(a^-aj  +  T-I  +  m  +  v  +  SA)  +  N  (2-60) 

A 

where 

^  =  carrier-phase  measurement  (cycles) 

A  =  carrier-phase  wavelength  (m) 

N  =  carrier-phase  integer  ambiguity  (cycles) 

The  remaining  terms  r,c,St^,St^^,T,I,m,v,SA  represent  the  same  error  types  as  in 

Equation  (2-59).  However,  the  magnitude  of  each  specific  error  type  in  a  code 
measurement  will  be  different  from  the  magnitude  of  the  same  error  type  in  the  carrier- 
phase  measurement.  In  addition,  the  ionospheric  term  is  the  same  magnitude  as  in  the 
pseudorange,  but  is  opposite  in  sign  because  the  ionosphere  delays  a  code  measurement 
but  advances  a  carrier-phase  measurement. 

A  GPS  receiver  cannot  distinguish  between  one  carrier  cycle  and  the  next.  The 
receiver  only  counts  the  change  in  the  number  of  cycles,  or  phase,  between  itself  and  the 
satellite.  The  initial  number  of  cycles  between  the  receiver  and  satellite  is  unknown,  or 
ambiguous.  To  make  the  carrier-phase  measurement  a  true  range  measurement,  an 
unknown  bias,  N ,  must  be  added  to  the  basic  receiver  measurement.  The  carrier-phase 
integer  ambiguity,  A ,  is  essentially  a  constant  unknown  integer-valued  bias  between  the 
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receiver  and  the  satellite,  and  is  different  for  each  receiver-satellite  measurement 


combination.  In  order  to  use  GPS  carrier-phase  measurements  for  positioning,  the 
carrier-phase  integer  ambiguity  must  be  estimated  or  resolved.  The  carrier-phase 
ambiguity  resolution  technique  used  in  this  research  is  discussed  in  Section  2.5.3. 

2.5.2.3  Single  Differencing 

In  addition  to  basic  pseudorange  and  carrier-phase  measurements.  Differential  GPS 
(DGPS)  techniques  use  measurements  created  by  differencing  combinations  of  receiver 
to  satellite  measurements.  The  concept  of  a  single-difference  measurement  is  shown  in 
Figure  2-2. 


Satellite  a 


Receiver  x  Receiver  y 

Figure  2-2.  Single  Differencing 

A  single-differenced  measurement  is  created  by  differencing  simultaneous 
measurements  between  two  receivers  and  one  satellite.  It  is  defined  as: 

^Pl=Pl-py  (2-61) 
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where  the  superscript  identifies  the  satellite  and  the  subscript  identifies  the  receiver  or 
receivers. 

PI  =  pseudorange  between  satellite  a  and  receiver  x 

Py  =  pseudorange  between  satellite  a  and  receiver  y 

Substituting  in  the  full  pseudorange  measurement  equation,  Equation  (2-59)  into  the 
above  equation  yields: 

-(r; -a;j+T; +i‘;+m‘;+v‘;+sA;) 

Rearranging  the  above  equation  to  combine  like  terms  yields 

Ap;  =  (r;  -  <) + c(a;.  -  -  a;,, ) + (t;  -  r;) 

+ (;;  -  /;) + {m-,  -».;)+ («;  -  «;> + (sa;  -  sa;) 

Because  the  measurements  differenced  in  Equation  (2-63)  are  from  the  same  satellite, 
the  satellite  clock  errors  and  the  SA  errors  are  eliminated,  (i.e.  a°  =0  and 


(2-62) 


(2-63) 


SA"  -  SA"  =  0 ).  Removing  these  terms  and  slightly  changing  notation  to  simplify  the 
expression  (allow  A  to  represent  a  single-difference)  the  final  single-differenced 
pseudorange  equation  becomes: 

Ap;  =  Ar;  +  +  AT"  +  A/^  +  +  Aw^  (2-64) 


The  single-differenced  pseudorange  measurement  is  really  a  measurement  of  the 
range  between  receivers  x  and  y  and  residual  errors.  The  residual  user  clock  error 
(Aca^^),  tropospheric  error  (AT^),  and  ionospheric  error  (AI^)  are  reduced  by  single 
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differencing,  while  the  multipath  errors  (Am^)  and  the  receiver  noise  errors  (Au^)  are 
slightly  increased  [37]. 

In  a  manner  similar  to  the  single-differenced  pseudorange  measurement,  the  single¬ 
difference  carrier-phase  measurement  can  be  written  as: 

(2-65) 


again  where  the  superscript  identifies  the  satellite  and  the  subscript  identifies  the  receiver 
and 

=  carrier-phase  measurement  between  satellite  a  and  receiver  x 

=  carrier-phase  measurement  between  satellite  a  and  receiver  y 

Substituting  in  the  full  carrier-phase  measurement  equation.  Equation  (2-60),  into  the 
above  equation  yields: 


AC 


-  + t; + n + < + c + ■501 + K 

-a;j+T;+i;  +/»;  +u; +s/i;)]+w; 


(2-66) 


Rearranging  the  equation  to  combine  like  terms  yields: 

AC  =|[(r;  -r;)+c«  -<»;.,)+(C -t;) 

A  (2-0/ ) 

+  (/;  -  /;)  +  K  -  <)  +  (^;  -  +  (*^4“  -  sa;)] + (n^  -  n;) 

Again,  the  satellite  clock  errors  and  SA  errors  are  eliminated: 


AC  =  jiK  +  Ac«;.  AT-;  -h  a;;  Am;  -f  AC )  aw;  (2-68) 


The  additional  term  in  Equation  (2-68),  AA^,  is  the  difference  in  the  carrier-phase 
ambiguities  between  the  two  receiver’s  measurements,  and  it  is  still  and  integer  value. 
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2.5.2.4  Differential  Corrections 

Code  based  DGPS  techniques  use  differential  corrections  to  improve  the  relative 
positioning  accuracy  of  GPS.  A  differential  correction  is  formed  by  subtracting  the 
calculated  range  to  a  satellite  from  the  measured  pseudorange: 

(2-69) 

where 

A/?“  =  the  differential  correction  between  receiver  x  and  satellite  a 

^  ^calc 

PI  =  the  measured  pseudorange  between  receiver  x  and  satellite  a 

=  the  calculated  range  between  receiver  x  and  satellite  a 

Differential  corrections  are  only  used  when  the  receiver’s  position  is  known,  as  it  is  in 
a  fixed  reference  receiver.  The  differential  correction  is  then  broadcast  or  used  in  post¬ 
processing  applications  to  correct  the  pseudorange  measurements  of  mobile  receivers. 
This  effectively  achieves  a  single-difference  measurement  between  the  user’s  receiver 
and  the  reference  receiver  as  shown  below: 

(2-70) 

p-  =r;+c(A  -AJ  +  Ar  +  AZ  +  Au  +  Am  (2-71) 

As  in  single-differencing,  the  satellite  clock  errors  and  SA  errors  are  eliminated  by  the 
differential  corrections. 
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2.5. 2. 5  Double  Differencing 


The  concept  of  double  differencing  is  shown  in  Figure  2-3  and  is  very  similar  to 
single  differencing.  A  double-difference  measurement  is  created  by  differencing  two 
single-differenced  measurements. 


Figure  2-3.  Double  Differencing 


The  double-differenced  measurement  is  defined  as: 

AV/>;  -  Ap;  -Ayo;  -Apj)  (2-72) 

Substituting  the  full  pseudorange  equation,  Equation  (2-59),  into  Equation  (2-72)  and 
combining  like  terms  yields: 

=(r;-r;)-(rj’-r;) 

+  (r;  -  r; )  -  )  (2-73) 

+  {ml  -mj) 

+  {SAi  -  sa;  )  -  (SA^^  -  SA'; ) 
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Again,  the  satellite  clock  errors  and  SA  errors  are  eliminated  as  well  as  the  receiver’s 
clock  errors.  Removing  these  terms  and  using  AV  to  represent  double  differencing  gives 


the  final  form  of  the  pseudorange  double-differencing  measurement  equation: 

=  AVr^"  +  AVr;"  +  AVZ  J  +  AVm  J  +  AV  v  J  (2-74) 

The  carrier-phase  double-difference  measurement  is  created  in  a  similar  manner 

-A«>*  =(Ar.  -Af*;)-(A«i,* -a«>;)  (2-75) 

Substituting  the  full  carrier-phase  measurement  equation,  Equation  (2-60),  into 
Equation  (2-75)  and  combining  like  terms  yields: 

Av<=|k--';)-(r/-r;) 

-t(r;-r;)-(r/-r;) 

+  (2-76) 

+  (m“  -mp-(m*  -m*) 

+  {SAl-SA';)-{SAl-SAl) 

Removing  the  eliminated  terms  and  using  the  double  difference  operator  gives  the 
final  form  of  the  double-differenced  carrier-phase  measurement  as: 

=  -  (AVr^  -I-  AVr^^  +  +  AVm J  +  AV (2-77) 

/I 

It  is  important  to  note  that  the  double-differenced  ambiguity  term  is  still  an 

integer. 
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2.5.2.6  Widelane  Measurements 


An  additional  carrier-phase  measurement  combination  is  the  widelane  measurement. 
A  widelane  measurement  is  a  linear  combination  of  LI  and  L2  phase  measurements 
defined  as: 


<l>wL=<l>a-<l>L2  (2-78) 

The  widelane  measurement  has  a  much  longer  wavelength  than  either  an  LI  or  L2 
measurement.  In  addition  to  its  long  wavelength,  the  widelane  measurement  remains  an 
integer  quantity.  The  longer  wavelength  of  widelane  measurements  decreases  the 
number  of  candidate  ambiguities  in  a  given  search  space,  making  the  widelane 
measurement  very  useful  in  ambiguity  resolution  routines.  Table  2-2  below  summarizes 
the  different  carrier-phase  measurements  combinations  and  associated  ambiguity  terms 
used  in  this  research. 


Table  2-2.  Phase  Measurement  Combinations 


Measurement 

Wavelength  (cm) 

Ambiguity 

LI 

19.03 

AVA^^j  (Integer) 

L2 

24.42 

AVNl2  (Integer) 

Widelane 

86.19 

(Integer) 

In  addition  to  widelane  measurements,  there  are  other  phase  measurement 
combinations  not  used  in  this  research,  but  useful  in  other  GPS  applications.  For  more 
detail,  the  reader  is  referred  to  references  [20,  37]. 
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2.5.3  Carrier-phase  Ambiguity  Resolution 

As  stated  earlier,  in  order  to  use  GPS  carrier-phase  measurements  for  positioning,  the 
integer  ambiguities  must  be  determined.  There  are  many  different  algorithms  used  to 
determine  integer  ambiguities.  However,  most  algorithms  perform  two  primary 
operations.  First,  they  determine  the  ambiguity  search  space  by  generating  sets  of 
ambiguities.  Second,  they  must  determine  which  set  of  those  generated  is  the  correct 
ambiguity  set. 

2.5.3.1  Ambiguity  Set  Generation 

The  ambiguity  set  generation  routine  used  in  this  research  is  a  combination  of  the 
LAMBDA  method,  developed  by  Teunissen  [40]  and  the  FASF  method  developed  by 
Chen  and  Lapachelle  [6]. 

Initial  floating-point  estimates  of  the  carrier-phase  ambiguities  and  their  associated 
covariance  matrix  are  generated  using  an  EKF.  However,  these  estimates  tend  to  be 
highly  correlated,  especially  for  short  time  periods,  causing  problems  in  determining  the 
correct  integer  ambiguity  set.  The  basic  idea  behind  the  LAMBDA  method  is  to 
decorrelate  the  highly  correlated  carrier-phase  ambiguities  by  applying  a  decorrelating 
ambiguity  transformation  known  as  the  “Z-transformation.”  Once  the  ambiguities  are 
decorrelated,  estimating  the  integer  carrier-phase  ambiguities  can  be  done  rapidly  and 
efficiently. 

The  Z-transform  used  in  the  LAMBDA  method  must  meet  certain  conditions  to  be 
effective  in  the  ambiguity  resolution  problem.  First,  to  avoid  resizing  the  search  space, 
the  Z-Transform  must  be  volume  preserving  (one-to-one  relation).  Second,  to  decorrelate 
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the  ambiguities,  it  must  reduce  the  product  of  ambiguity  variances.  Finally,  to  preserve 
the  integer  nature  of  ambiguities,  the  Z-Transform  must  have  integer  elements.  The  Z- 
Transform  used  in  the  LAMBDA  method  is  defined  by: 

z  =  Z-x  ,  z  =  Z-x  ,  Q,  =ZQiZ’’  (2-79) 

where 


z,  X  =  transformed  and  untransformed  ambiguities 

Z  =  Z-Transformation  matrix 

z,x  =  transformed  and  untransformed  ambiguity  estimates  ambiguity  estimates 

Qz  »Qx  =  transformed  and  untransformed  ambiguity  covariance  matrices 

The  Z-Transform  (for  a  two  dimensional  case)  is  constructed  by  the  following 
equations.  Assuming 


then 


(2-80) 


1 

and  Z,  = 

1 

O' 

0  1 

_int(-<T2it7“^) 

1 

(2-81) 

where  either  the  upper  diagonal  form  Zj,  or  the  lower  diagonal  form  Z2,  of  the 
transform  may  be  used.  In  the  above  equation,  int  refers  to  rounding  to  the  closest  integer 
value,  which  is  necessary  to  preserve  the  integer  nature  of  the  carrier-phase  ambiguities. 
Without  rounding  the  elements  to  integer  values.  Equation  (2-81)  would  fully  decorrelate 
the  ambiguities,  but  the  ambiguities  would  no  longer  be  integer  values. 
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An  efficient  method  to  calculate  higher  order  Z-transforms,  used  in  this  thesis,  was 
developed  by  Rizos  and  Han  [38].  The  first  step  in  the  method  is  to  perform  an  upper 
triangular  factorization  of  the  ambiguity  estimate  covariance  matrix  by: 

Q,=U,D„U[  (2-82) 

Next,  the  integer  matrix  is  calculated  by: 

=[int(U,)r  (2-83) 

Then  an  intermediate  covariance  matrix  is  calculated  by: 

(2-84) 

The  second  step  is  to  perform  a  lower  triangular  factorization  of  the  intermediate 
covariance  matrix  by: 

Qz,  (2-85) 

Next,  the  integer  matrix  is  calculated  by: 

=[int(L,)]-‘  (2-86) 

Then  an  intermediate  covariance  matrix  is  calculated  by: 

(2-87) 

This  procedure  is  iterated  until  both  integer  matrices  [int(Uj)]  and  [int(Lj)]  become 
identity  matrices.  The  overall  Z-transformation  is  obtained  by: 

Z  =  Z,^_Z,^_,...Z^Z,^Z^Z,,  (2-88) 

The  Z-Transform  calculated  in  (2-88)  is  used  to  transform,  or  decorrelate,  the 
floating-point  ambiguity  estimates  and  the  associated  covariance  matrix  before 
determining  the  ambiguity  search  space.  Actual  generation  of  candidate  ambiguity  sets 


for  this  thesis  is  accomplished  using  the  transformed  ambiguities  via  the  FASF  method, 
but  can  be  done  using  many  of  the  routines  described  in  Section  1. 

The  FASF  method  starts  by  determining  the  search  space  for  the  best-known  single 
ambiguity.  The  search  space  is  determined  by  using  a  predefined  multiple  (k)  of  the 
standard  deviation  of  the  ambiguity.  The  search  space  for  the  ambiguity  is: 

x.-ka.<WiM^<x,+ko,  (2-89) 

where 

Jc„  =  floating-point  estimate  of  the  n^^  ambiguity 

th 

=  Standard  deviation  of  the  n  ambiguity 

For  each  possible  integer  value  in  this  search  space,  the  integer  ambiguity  is  assumed 
to  be  the  correct  value.  Then  new  conditional  covariance  and  ambiguity  estimates  are 
calculated,  conditioned  on  the  fact  that  the  first  ambiguity  is  known.  This  process  shrinks 
the  covariance  terms  of  the  remaining  ambiguities.  Next,  the  best-known  integer  value 
from  the  remaining  ambiguities  is  selected,  and  the  process  repeated.  The  method  is 
recursive  and  will  either  result  in  a  valid  ambiguity  set  or,  if  the  conditional  covariance 
doesn’t  allow  for  a  valid  integer  somewhere,  the  process  will  result  in  an  invalid  set  that 
will  be  removed  from  the  candidate  ambiguity  sets.  The  new  conditional  state  estimate 
and  covariance  are  calculated  by: 

x  =  x-p„(x„-S/AN,J/ctI  (2-90) 

P.=P.-(pX)/^„^  (2-91) 
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where 


X  =  unconditioned  estimated  parameter  vector 
X  =  estimated  parameter  vector  conditioned  upon 

P-  =  covariance  matrix  conditioned  upon  =  VAA^i„, 
p„  =  n'*  column  of  the  covariance  matrix  P- 

cr^  =  scalar  variance  of  the  parameter  (taken  from  the  diagonal  of  P^ ) 

The  FASF  method  can  be  used  to  generate  ambiguity  sets  in  either  the  normal  or  Z- 
transformed  search  space.  However,  the  candidate  ambiguity  sets  must  be  transformed 
back  to  the  original  search  space  before  being  used  in  GPS  measurements. 


2.5.3. 2  Ambiguity  Set  Determination 

Most  carrier-phase  ambiguity  resolution  algorithms  monitor  the  measurement 
residuals  for  each  candidate  ambiguity  set,  by  statistically  comparing  possible  solutions, 
until  the  correct  ambiguity  set  becomes  evident.  One  method  of  comparing  the  residuals 
is  by  the  ratio  test: 


ratio  = 


r^R~V(2"‘^best) 

r^R“V(best) 


(2-92) 


where 

r  =  measurement  residuals 
R  =  measurement  covariance 

The  ambiguity  set  is  determined  to  be  the  correct  ambiguity  if  the  ratio  is  consistently 
greater  than  a  given  value,  typically  2. 
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This  is  a  common  method  of  detemiining  the  correct  ambiguity,  however  it  is  not  the 
method  used  in  this  research.  The  method  used  in  this  research  is  to  base  each  one  of  the 
elemental  filters  in  an  MMAE  on  a  different  ambiguity  set.  The  conditional  probability 
of  each  filter  being  based  on  the  correct  ambiguity  set  is  then  calculated  by  Equation 
(2-55),  with  the  modifications  outlined  in  Section  3.5.2.  Based  on  its  convergence 
properties,  the  MMAE  will  converge  to  the  correct  ambiguity  set,  provided  the  correct 
ambiguity  set  is  used  in  one  of  the  elemental  filters.  If  the  correct  ambiguity  set  is  not 
available  in  one  of  the  elemental  filters  the  MMAE  will  converge  to  the  elemental  filter 
with  the  ambiguity  set  “closest”  to  the  correct  ambiguity  set  in  the  Barum’s  distance- 
measure  sense  [4]. 

2.6  Summary 

This  chapter  presented  the  introduction  and  theory  of  linear  and  Extended  Kalman 
filters,  the  MMAE  and  GPS  (including  carrier-phase  ambiguity  resolution  techniques). 
Chapter  3  will  provide  detailed  descriptions  of  the  models  and  algorithm  development 
investigated  in  this  thesis. 
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3  Methodology  and  Algorithm  Development 


3.1  Overview 

This  chapter  provides  a  complete  description  of  the  overall  algorithm,  system  models, 
and  techniques  used  in  this  thesis.  First,  an  overall  algorithm  description  is  presented, 
followed  by  a  description  of  each  main  subroutine.  Detailed  descriptions  of  the  system 
states,  models  and  measurements  equations  for  the  Kalman  filters  used  in  this  research 
are  then  presented.  Finally,  specific  techniques  are  presented  that  are  key  to  the 
understanding  of  this  research,  including  changes  made  to  the  theory  of  Chapter  2. 

3.2  Overall  Algorithm  Structure 

The  overall  structure  of  the  algorithm  developed  in  this  thesis  is  shown  in  Figure  3-1. 
The  algorithm  can  easily  be  divided  into  three  main  areas:  a  floating-point  differential 
filter,  carrier-phase  ambiguity  set  generation  (Z-transform,  FASF  and  inverse  Z- 
transform),  and  the  MMAE.  Each  one  of  these  areas  is  fully  described  in  the  following 
sections. 
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Figure  3-1.  Overall  Algorithm  Description 

3.3  Floating-Point  Differential  GPS  Kalman  Filter 

The  initial  step  in  the  post-processing  algorithm  is  a  floating-point  DGPS  Kalman 
filter,  and  it  follows  the  general  development  of  the  filter  found  in  reference  [36],  This 
filter  is  a  DGPS  navigation  EKF  designed  to  calculate  the  user’s  position,  velocity. 


_ I _ I _ 

Conditional  Probability 
Calculation 


acceleration,  and  carrier-phase  ambiguity  estimates.  The  filter  uses  both  differentially 


corrected  code  measurements  and  double-differenced  carrier-phase  measurements.  The 
main  purpose  of  this  initial  filter  is  to  generate  the  carrier-phase  ambiguity  estimates  and 
covariances  used  by  the  ambiguity  set  generation  routine.  However,  the  output  of  this 
filter  is  also  used,  as  necessary,  to  initialize  or  reinitialize  the  elemental  filters  in  the 
MMAE. 

A  few  aspects  are  critical  to  implementing  this  DGPS  filter.  First,  it  is  assumed  that 
the  remote  (or  airborne)  receiver  has  access  to  the  same  satellites  as  the  reference 
receiver.  Second,  it  is  also  assumed  for  modeling  purposes  that  the  differential 
corrections  are  all  timed  correctly  to  correspond  to  the  correct  remote  receiver 
pseudorange  measurement.  Last,  when  processing  DGPS  measurements,  it  is  assumed 
that  the  differential  corrections  have  been  applied  to  the  raw  pseudorange  measurements 
from  the  remote  receiver. 

All  of  the  above  assumptions  are  easily  and  properly  invoked  because  the  algorithm 
developed  in  this  thesis  is  a  post-processing  algorithm.  A  simple  data  preprocessing 
routine  was  used  to  coordinate  visible  satellites  and  measurement  times  between  the 
reference  and  remote  receivers  prior  to  the  data  being  used  in  the  floating-point  DGPS 
Kalman  filter.  In  addition,  the  data  preprocessing  routine  was  used  to  detect  cycle  slips, 
determine  a  base  satellite,  and  to  apply  tropospheric  corrections  to  raw  pseudorange 
measurements.  The  data  preprocessing  routine  does  not  in  any  way  affect  the  functioning 
of  the  overall  algorithm.  It  simply  organizes  the  available  data  to  improve  the  efficiency 
of  the  floating-point  DGPS  Kalman  filter  and  the  MMAE. 
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3.3.1  Differential  GPS  Model  Equations 


The  floating-point  Differential  GPS  Kalman  filter  is  an  ll-h(n-l)  state  filter  (n  is  the 

number  of  satellites  in  view)  with  a  linear  dynamics  model.  The  filter’s  first  nine  states 

consist  of  three  position,  three  velocity  and  three  acceleration  states,  based  on  a  First 

Order  Gauss  Markov  Acceleration  (FOGMA)  model. 

The  position  and  velocity  states  dynamics  are  defined  as: 

x^  =  X  =x^  x^=Y  =  x^  Xg  =  Z  =  Xg 
x^  —  it  =  X.J  Xj— T  —  Xg— Z  —  JCg 

To  handle  unknown  accelerations,  the  acceleration  states  are  modeled  as  First  Order 
Gauss  Markov  processes.  The  acceleration  states  dynamics  are  defined  as: 

^7  =(-V^aK+w«,(0 

i:,=(-l/r,K+>v„^(0  (3-1) 

The  acceleration  states  process  noises  are  defined  as: 


if  +  T)}=  +  r)}=  (t  +  T)} 


^^  =  q^S{T) 


The  time  constant,  ,  and  the  mean  squared  value,  erf ,  are  user-defined  tuning  values, 

set  according  to  the  expected  accelerations  and  time  correlations. 

The  largest  error  sources  in  a  GPS  measurement  are  due  to  clock  errors  in  the  user’s 
receiver.  Therefore,  two  additional  states  were  added  to  the  first  nine  states,  to  represent 
the  user’s  clock  bias  and  clock  bias  drift.  The  clock  bias  state  and  drift  state  dynamics 


are  defined  as: 


(3-3) 


X,o 

•*11  = 


The  clock  bias  and  drift  states  process  noises  are  defined  as: 


E{v^ik,  (t>dk,  (t  +  ^)}=  (t  +  ^■)}=  (3-4) 


The  floating-point  DGPS  filter  uses  both  single-differenced  (differentially  corrected) 
code  measurements  and  double-differenced  carrier-phase  measurements.  The  double- 
differenced  carrier-phase  measurements  may  be  LI,  L2,  or  widelane  measurements. 
Note,  widelane  measurements  are  used  for  the  majority  of  this  research.  In  order  to  use 
double-differenced  carrier-phase  measurements,  the  filter  must  also  estimate  the  carrier- 
phase  ambiguities.  By  definition,  the  double-difference  carrier-phase  integer  ambiguities 
are  constant  biases  in  the  carrier-phase  measurements.  However,  due  to  errors  in  the  GPS 
measurements,  if  these  states  were  modeled  as  constant  biases  the  filter  may  converge  to 
the  wrong  values,  with  no  way  to  correct  the  bias.  A  bias  in  the  ambiguity  estimates  will 
result  in  a  constant  bias  in  the  position  estimate  of  the  overall  algorithm.  Therefore,  these 
states  were  modeled  as  random  walks,  and  an  additional  n-\  states  were  added  to  the 
filter  to  represent  the  double-differenced  carrier-phase  ambiguities.  The  double- 
differenced  carrier-phase  ambiguity  states  dynamics  are  defined  as: 


(3-5) 


'-12  - 


with  process  noises  defined  as: 


•  •  =  (^)^VAW” 


(3-6) 
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The  overall  state  vector  for  the  floating-point  filter  is  defined  by: 
x=lx  Y  Z  X  Y  Z  X  Y  Z  cdt^^  VAiV*^  •••  (3-7) 

where 

x^-X=  ECEF  X  position  (m) 

^2  =  Y  =  ECEF  Y  position  (m) 

X3  =  Z  =  ECEF  Z  position  (m) 

X4  =  X  =  ECEF  X  velocity  (m/s) 

X5  =  Y  =  ECEF  Y  velocity  (m/s) 
jCg  =  Z  =  ECEF  Z  velocity  (m/s) 

X7  =  X  =  ECEF  X  acceleration  (m/s^) 

Xg  =  Y  =  ECEF  Y  acceleration  (m/s^) 

JC9  =  Z  =  ECEF  Z  acceleration  (m/s^) 

jCjo  =  cS  =  user’s  (remote)  and  reference  clock  bias  (m) 

x,j  =cS  i^ii^  =  user’s  (remote)  and  reference  clock  drift  (m/s) 

Xj2  =  VA/V*^  =  double-difference  phase  ambiguity  between  satellites  1  and  2  (cycles) 

Xj3  =  =  double-difference  phase  ambiguity  between  satellites  1  and  3  (cycles) 

^ii+(n-i)  “  VAAT'"  =  double-difference  phase  ambiguity  between  satellites  1  and  n  (cycles) 

Because  GPS  measurements  are  in  the  Earth  Centered  Earth  Fixed  (ECEF)  reference 
frame,  all  of  the  calculations  in  the  floating-point  filter  (and  the  MMAE)  were  conducted 
in  the  ECEF  frame. 
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The  linear  state-space  differential  equation  for  the  floating-point  filter  in  the  form: 

x(t)=F(t)x(t)+ G(t)w(t)  (3-8) 


is  fully  represented  by: 
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In  Equation  (3-9),  Ta  is  the  time  constant  associated  with  the  FOGMA  acceleration 
states.  For  this  system,  the  matrix  G(t)  is  an  identity  matrix.  Although  the  floating-point 
filter  does  not  model  all  of  the  possible  error  sources  in  a  DGPS  application,  enough  are 
modeled  to  be  accurate  enough  for  the  purposes  of  this  research. 
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The  dynamics  driving  noise  matrix  Q  defined  by: 


(t  +  t)}=  Q(t)j(T)  (3-10) 


must  also  be  specified  for  the  filter.  The  matrix  Q  used  in  this  research  was  given  by: 
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The  floating-point  filter  time  constant,  mean  squared  value,  and  noise  values  used  in  this 
thesis,  for  both  ground  and  flight  test,  are  specified  below  in  Table  3-1. 


Table  3-1.  Floating-Point  Filter  Dynamics  Driving  Noise  Values 


Term 

Deflnition 

Value 

Mean  squared  value 

(12.25  m/sec')" 

Ta 

Acceleration  time  constant 

3  sec 

Acceleration  noise 

100  m^  /sec^ 

^clki 

Clock  bias  noise 

0.036  m^/ sec 

^Clk2 

Clock  Drift  noise 

0.141mVsec' 

Ambiguity  noise 

1 . 1 X 10^  cycles^/sec 
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The  initial  conditions  for  the  position  states  {x^,X2,x^)  and  the  clock  bias  state 
are  determined  using  a  single  point  GPS  positioning  subroutine.  The  double-differenced 
ambiguity  states  are  initialized  using  code-carrier  difference  calculations, 

by  the  following  equation: 

x^„  =  VAN^"  »  YAf"  -  (3-12) 

A 

All  remaining  states  were  initialized  to  zero,  making  the  initial  state  vector: 
x(to)  =  [Xo  To  ^0  0  0  0  0  0  0  0  VAAq"  VAAf  •••  (^-13) 

The  initial  covariance  matrix  for  the  filter  is  assumed  to  be  a  diagonal  matrix  with 
elements  shown  below.  All  of  the  initial  cross-covariance  terms  set  to  zero.  The  initial 


covariance  matrix 

for  the  filter 

is  given  by: 
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The  initial  covariance  values  used  in  this  thesis  are  specified  below  in  Table  3-2. 


Table  3-2.  Floating-Point  Filter  Initial  Covariance  Values 


Term 

Deflnition 

Value 

(7^ 

Position  state  variance 

(100m)' 

2 

Velocity  state  variance 

(400m /^)' 

<7^ 

X^,X^X^ 

Acceleration  state  variance 

(20  m/ 5')' 

■*10 

Clock  bias  variance 

(lOOm/s^y 

Xii 

Clock  drift  variance 

(1000  m/s^y 

^2  _2 
•*^12  Xi„ 

Ambiguity  variance 

i^cyclesf 

3.3.2  Differential  GPS  Measurement  Model 

The  floating-point  DGPS  Kalman  filter  uses  a  nonlinear  measurement  model 
consisting  of  differentially  corrected  code  measurements  and  double-differenced  carrier- 
phase  measurements  resulting  in  a  2n-l  (n  is  the  number  of  satellites  in  view) 
measurement  vector  specified  by: 

pI„  -  Pl,r  (3-14) 


where 

Plorr  -  differentially  corrected  code  measurement  between  the  receiver  and  satellite  1 
Plorr  -  differentially  corrected  code  measurement  between  the  receiver  and  satellite  2 

=  differentially  corrected  code  measurement  between  the  receiver  and  satellite  n 
VA(z!>^  =  double-differenced  carrier-phase  measurement  satellites  1  and  2 
VA(^^  =  double-differenced  carrier-phase  measurement  satellites  1  and  3 

VA^^]^  =  double-differenced  carrier-phase  measurement  satellites  1  and  n 
To  be  used  in  an  extended  Kalman  filter  the  measurement  model  must  be  in  the  form: 
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z{t,)=h[x(ti),t,]+  v(f,.) 


(3-15) 


From  Section  2.5.2.4,  differentially  corrected  code  measurements  are  defined  as: 

=r°  +c(a-St^  +  AT  +  M  +  Av  +  ^m  (3-16) 

By  combining  terms,  the  above  equation  becomes: 

PL  +  (3-17) 

where  cSt  =  c{6ty  -  6t^)  and  =  AF  -I-  A/  -t-  Au  +  Am .  In  Equation  (3-17),  the  noise 

term,  ,  is  modeled  as  a  white  noise.  Even  though  the  individual  residual  differential 

error  terms  (AT ,AI ,Av  ,Am)  are  not  white  noises,  we  assume  the  combination  of  the 
four  terms  can  be  sufficiently  modeled  as  a  white  noise  due  to  the  central  limit  theorem. 
Expanding  the  range  term,  to  incorporate  satellite  position  and  receiver  position,  and 
substituting  state  variables  for  the  appropriate  terms  gives  the  following  equation: 

+  x.o-h  Vp  (3-18) 

This  measurement  equation  must  be  linearized  before  being  used  in  the  gain  and 
covariance  calculations  of  an  EKF.  The  partial  derivative  matrix  H  is  defined  as: 
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The  individual  partial  derivatives  for  one  row  of  the  differentially  corrected  code 


measurements  are: 


Sx^ 


Xsv-Xi 


x=Mtn  )^  +  (>'.v  -  ^2  )^  +  (^^  -  ^3  F 


(3-20) 


Sh[x,ti] 


Si[x,tj] 


ysv-X2 


<=Ht7)  t.v  -^l)^  +(}'.V  -^2)^  +U.V  --^3)^1'^ 


&C. 


2.v-^3 


x=x(ti  ) 


[(^5V  -^1)^  +(j^v  -^2)^  +(^JV 


=  ^2  (3-21) 


=  (3-22) 


Sh[X,  t;  ] 


&c, 


10 


=  1 


(3-23) 


x=Hh) 


These  can  be  combined  into  a  vector  representing  one  row  of  the  H  matrix  as: 

H‘'  =  le“  0  0  0  0  0  0  1  0  Oj  (3-24) 


where 


e“  =  [^1  62  ^3]  =  unit  line  of  sight  vector  between  the  receiver  and  satellite  a. 


The  double-differenced  carrier-phase  measurements  are  used  in  a  similar  manner. 
From  Section  2.5.2.5,  the  double-differenced  carrier-phase  measurements  are  defined  by: 

=  -  (AVr^"  +  ^S/T^  +  AV/ +  AVm J  -h  AV )  +  AVAT J  (3-25) 

A 

The  A,  term  in  Equation  (3-25),  represents  the  wavelength  of  the  carrier-phase 
measurement  and  will  be  different  depending  on  the  type  of  double-differenced  carrier- 
phase  measurement  (LI,  L2,  Widelane)  being  used  in  the  filter. 
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By  expanding  the  double-differenced  range  term  and  combining  the  residual 
differential  error  terms  the  above  equation  becomes: 

AV<=j[r;-r;-(r;-r‘)]+  (3-26) 

In  Equation  (3-26),  is  modeled  as  a  white  noise,  and  is  the  combination  of  the 
residual  error  terms  in  the  double-differenced  equation,  i.e. 

V =  —  .  Again,  even  though  the  individual 

A 

residual  error  terms  (AVr^,AV/^*,AV/n^,AVv^)  are  not  white  noises,  as  explained 

previously,  the  combination  of  the  four  terms  can  be  sufficiently  modeled  as  a  white 
noise.  By  expanding  the  range  terms  and  substituting  state  variables  for  the  appropriate 
terms  gives  the  following  equation; 

=|[(a:" -x,)^  +(y‘ -x,)‘+(z‘  -x,)^]'" 

-kx''  -x,f  +(y‘'  -x^f  -hiz"  -x,y]'^  (3-27) 

+  J  -'•;}+ A 

The  individual  partial  derivatives  for  one  row  of  the  double-differenced  widelane 
carrier-phase  measurements  are: 
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Sh[x,ti]  _  1  j _ x‘‘  -x^ _ ^ 

-x,f+(y‘ 

_j_| _ x'^  -X^ _ 1  (3-28) 

-l  1 -  ;c, )"  +  (y*  -  X,)'  +  (z*  - r 

4k -d 


(3-29) 


where 

e"  =  [gj  6^]=  unit  line  of  sight  vector  between  the  receiver  and  satellite  a. 
e*  =  [(2j  ^2  ^3]-  line  of  sight  vector  between  the  receiver  and  satellite  b. 
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The  additional  “1”  in  the  above  equation  is  placed  in  the  column  corresponding  to  the 
appropriate  ambiguity  state. 

The  entire  measurement  matrix  H  is  an  (1  l+(n-l))  by  (2n-l)  matrix  defined  by; 
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The  measurement  error  covariance  matrix  R  defined  by: 


R(r,.)  for  t,=tj 
0  for  t,.  tj 


(3-33) 


(3-34) 


must  also  be  specified  for  the  filter.  The  matrix  must  consider  code  measurement  noise 
variances,  phase  measurement  noise  variances,  the  code  measurement  noise  covariances, 
the  phase  measurement  noise  covariances,  and  the  cross-covariance  of  code  measurement 
and  phase  measurement  noise.  For  this  research,  the  code  measurement  noise  covariance 
was  assumed  to  be  zero  based  on  the  assumption  that  code  measurement  errors  are 
uncorellated  between  satellites,  and  the  baseline  distances  are  short  enough  to  cancel  the 
majority  of  tropospheric  errors.  In  addition,  the  cross-covariance  of  code  and  phase 
measurement  noises  was  assumed  to  be  zero  based  on  the  assumption  that  double- 
differenced  tropospheric  errors  are  negligible  and  that  phase  multipath  and  noise  errors 
are  uncorrelated  between  measurements.  Values  for  the  other  measurement  error  terms 
are  listed  in  Table  3-3.  The  R  matrix  is  a  (2n-l)-by-(2n-l)  matrix  specified  by: 
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where 

r^i  =  differentially  corrected  code  measurement  variance 

VA<»»  “  double-differenced  carrier-phase  measurement  variance 
r^A^.j  vA^*  “  double-differenced  carrier-phase  measurement  covariance 

The  floating-point  filter  measurement  error  values  used  in  this  research  are  listed  below 
in  Table  3-3. 


Table  3-3.  Floating-Point  Filter  Measurement  Covariance  Error  Values 


Term  I  DeHnition  I  Value 


Differentially  corrected  code  variance  error 


'vA«>ff,vA(»'J  Double-differenced  carrier-phase  variance  error  0.004m^ 
^  ,vAf‘  Double-differenced  carrier-phase  covariance  error  0.002m^ 


Value 


2.6m^ 


0.004m^ 


VA0'-'  ,VA^‘ 
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3.3,3  Discrete-Time  Models 


The  equations  derived  above  must  be  converted  into  difference  equations  in  order  to 
be  implemented  in  a  Kalman  filter.  The  filter  model  equations  must  be  in  the  form; 

(3-36) 

where 

£{w,}=0 

(3-37) 

Mw,(r^.)w^(r^)}=0, 

The  state  transition  matrix,  ^(At) ,  is  calculated  by: 

=  =  (3-38) 

where  At  =  -4 ,  resulting  in  the  matrix: 

'lOOAtO  OAOOOOOOOOl 
OlOOAtO  OAOOOOO 
0010  OArOOAOOOO 
0001  0  OBOOOOOO 

0000  1  OOBOOOOO 

0000  0  1  OOBOOOO 

0000  0  ocoooooo 

d)(At)  =0000  0  0  ocooooo 

0000  0  0  oocoooo 

0000  0  0  0001  At  00 

0000  0  0  0000100 

0000  0  00000010 

0000  0  0  0000001 


[0  000  0  0  000000001 

where 


•  0 
•  0 
•  0 
•  0 
•  0 

•  0  (3-39) 

•  0 
•  0 
•  0 
•  0 
•  0 
.  0 
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The  discrete  dynamics  driving  noise  matrix  is  calculated  by: 
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3.4  Carrier-phase  Ambiguity  Set  Generation 

The  carrier-phase  ambiguity  set  generation  routine  consists  of  the  Z-transform,  FASF 
and  inverse  Z-transforms  subroutines.  The  double-differenced  ambiguity  states  and  their 
associated  covariance  output  from  the  floating-point  filter  are  used  in  the  carrier-phase 
ambiguity  set  generation  routine. 

The  first  step  in  generating  candidate  ambiguity  sets  is  to  generate  the  Z-transform. 
The  Z-transform  is  the  linear,  volume  preserving,  integer  transform  described  in  Section 
2.5.3. 1.  The  Z-transform  used  in  this  algorithm  is  based  on  the  higher  order 
transformations  detailed  in  Equations  (2-82)  through  (2-88).  Once  generated,  the  Z- 
transform  is  used  to  decorrelate  the  floating-point  DGPS  filter  ambiguity  estimates  before 
being  used  in  the  FASF  subroutine. 

These  transformed  estimates  are  then  used  by  the  FASF  subroutine  also  described  in 
Section  2.5.3. 1  and  Equations  (2-89)  through  (2-91).  The  FASF  subroutine  generates  a 
user-defined  number  of  candidate  ambiguity  sets,  picked  to  be  250  sets  for  this  thesis. 
Then  the  routine  picks  the  K  best  candidate  sets  (where  K  is  the  number  of  elemental 
filters  used  in  the  MMAE,  between  15  and  100,  as  described  in  Section  3.5)  by  using  the 
ratio  test  in  Equation  (2-92).  The  K  candidate  ambiguity  sets  are  transformed  back  into 
the  original  measurement  space  using  the  inverse  of  the  Z-transform  before  being  used  in 
the  elemental  filters  of  the  MMAE. 

New  candidate  ambiguity  sets  are  only  generated  and  provided  to  the  elemental  filters 
in  the  MMAE  under  certain  conditions: 
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□  During  initialization,  candidate  ambiguity  sets  are  generated  after  a  single  propagate 
and  update  cycle  of  the  floating-point  filter. 

□  New  candidate  ambiguity  sets  are  generated  after  uncorrected  cycle  slips  in  either  the 
reference  or  remote  receivers. 

□  New  candidate  ambiguity  sets  are  generated  when  a  satellite  is  lost  from  view  or  a 
new  satellite  comes  into  view. 

The  above  conditions,  although  simplistic  and  overly  restrictive,  were  implemented  to 
validate  the  “proof  of  concept”  algorithm  developed  in  this  thesis.  An  operational  system 
would  have  advanced  logic  routines  capable  of  more  efficiently  handling  cycle  slips  and 
satellites  changes. 

3.5  Multiple  Model  Adaptive  Estimator  Design 

As  stressed  throughout  this  thesis,  the  MMAE  design  is  the  key  to  determining  the 
correct  carrier-phase  ambiguity  set,  and  the  generation  of  accurate  GPS  navigation  data. 
The  MMAE  design  used  in  this  thesis  is  the  design  described  in  Section  2.4,  with  minor 
changes  to  the  conditional  probability  generation  calculations.  In  addition,  a  filter 
pruning  routine  was  added  to  allow  the  MMAE  to  converge  to  the  true  fixed-point 
solution. 
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3.5.1  Elemental  Filter  Design 

The  elemental  filters  in  the  MMAE  are  nine-state  (3  position,  3  velocity,  and  3 
acceleration)  FOGMA  filters,  described  below.  The  additional  clock  or  ambiguity  states 
used  in  the  floating-point  filter  are  not  necessary  because  the  elemental  filters  only  use 
double-difference  carrier-phase  measurements.  All  of  the  elemental  filters  are  identical, 
i.e.  they  use  the  same  system  dynamics,  the  same  states,  the  same  process  noises  and  the 
same  tuning  values.  This  means  there  is  only  one  dynamics  “model”  in  the  Multiple 
Model  Adaptive  Estimator  (MMAE).  The  only  difference  between  elemental  filters  is 
that  each  filter  “corrects”  the  double-differenced  carrier-phase  measurements  based  on  a 
different  double-differenced  ambiguity  set.  At  every  sample  time,  each  elemental  filter 
receives  (n-1)  double-difference  carrier-phase  measurements,  corrects  them  according  to 
its  own  double-differenced  ambiguity  set.  Each  elemental  filter  then  uses  the  “ambiguity 
corrected”  double-difference  carrier-phase  measurements  to  generate  state  estimates  and 
residuals. 

Because  each  elemental  filter  uses  a  fixed  ambiguity  set  to  correct  the  carrier-phase 
measurement,  at  every  sample  time  each  elemental  filter  generates  a  fixed-point  solution. 
The  elemental  filter’s  state  estimates  and  residuals  are  used  in  the  conditional  probability 
calculations  to  determine  which  elemental  filter’s  solution  is  the  true  fixed-point  solution. 

3.5.1.1  Elemental  Filter  Model  Equations 

The  state  vector  for  the  elemental  filters  are  given  by: 

x  =  [x  Y  Z  X  Y  Z  X  Y  z]  (3-42) 

where 
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jc,  =  X  =  ECEF  X  position  (m) 

^2  =  Y  =  ECEF  Y  position  (m) 

Xj  =  Z  =  ECEF  Z  position  (m) 

X4  =  X  =  ECEF  X  velocity  (m/s) 

X5  =  Y  =  ECEF  Y  velocity  (m/s) 

Xg  =  Z  =  ECEF  Z  velocity  (m/s) 

X2  =  X  =  ECEF  X  acceleration  (m/s^) 

Xg  =  Y  =  ECEF  Y  acceleration  (m/s^) 

X9  =  Z  =  ECEF  Z  acceleration  (m/s^) 

The  dynamics  and  process  noises  for  the  elemental  filter  were  derived  in  the  same 
manners  as  the  first  nine  states  of  the  floating-point  filter.  The  linear  state  space 


differential  equation  for  the  elemental  filters  is  fully  represented  by: 
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Again,  Ta  is  the  time  constant  associated  with  the  FOGMA  acceleration  states.  For 
this  system,  the  matrix  G(t)  is  an  identity  matrix.  The  dynamics  driving  noise  vector 
matrix  Q  must  also  be  specified  for  the  elemental  filters  and  is  given  by: 
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The  elemental  filter  time  constant,  mean  squared  value,  and  noise  values  used  in  this 
thesis,  for  both  ground  and  flight  test  are  specified  below  in  Table  3-4.  The  mean 
squared  value  and  hence  the  acceleration  noises  were  increased  when  processing  flight 
test  data  to  account  for  unknown  accelerations  not  present  during  the  static  ground  tests. 


Table  3-4.  Elemental  Filter  Dynamics  Driving  Noise  Values 


Term 

DeHnition 

Ground  Test  Value 

Flight  Test  Value 

Mean  squared  value 

(3.873m/sec")" 

(12.25m/sec^)^ 

T. 

Acceleration  time  constant 

3  sec 

3  sec 

Qa 

Acceleration  noise 

lOm^/sec^ 

100  m^/sec^ 

At  the  start  of  the  algorithm,  or  whenever  the  MMAE  is  reset,  the  elemental  filters  are 

* 

initialized  to  the  output  of  the  floating-point  filter. 


(3-45) 

(3-46) 
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3.5. 1.2  Elemental  Filter  Measurement  Model 


The  elemental  filters  use  a  discrete-time,  nonlinear  measurement  model  consisting 
only  of  double-differenced  carrier-phase  measurements.  Before  being  used  by  the 
elemental  filters,  the  double-differenced  carrier-phase  measurements  are  corrected  by  the 
candidate  ambiguity  set  being  tested  by  the  elemental  filter.  This  results  in  an  n-\ 
measurement  vector  specified  by: 


i{ti) 


(3-47) 


where 

(VA^^  )  =  ambiguity  corrected  double-differenced  carrier-phase 
measurement  between  satellites  1  and  2 

(VA^^  )  =  ambiguity  corrected  double-differenced  carrier-phase 
measurement  between  satellites  1  and  3 

(vA^^  )  =  ambiguity  corrected  double-differenced  carrier-phase 
measurement  between  satellites  1  and  n 


The  ambiguity-corrected  double-difference  carrier-phase  measurements  are  given  by: 

-  VAiV^J  (3-48) 


where 


VA^^*"  =  the  normal  (or  uncorrected)  double-differenced  carrier-phase 
measurement  between  satellites  a  and  b 
WAN^  =  the  double-differenced  candidate  ambiguity  term  between 
satellites  a  and  b 

The  candidate  double-differenced  ambiguity  set  (n-1  ambiguities)  being  tested  by  an 
elemental  filter  is  made  of  candidate  ambiguity  terms,  W AN ^ ,  one  ambiguity  for  each  of 
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the  double-differenced  measurements.  The  candidate  ambiguity  terms,  ,  are  by 

definition  constant  values,  but  are  different  constant  values  for  each  elemental  filter.  In 
fact,  the  difference  in  the  candidate  ambiguity  sets  is  the  only  difference  between  the 
elemental  filters. 

By  substituting  in  the  full  double-differenced  widelane  erarrier-phase  measurement 
equation  into  Equation  (3-48),  the  ambiguity  corrected  measurement  becomes: 

=  -(AVr^'’  -h  AVr^"  -H  Avrj  +  +  AVv^)  +  AWN^  -  AVN^^  (3-49) 

A 

As  outlined  in  Section  3.3.2,  the  measurement  model  must  be  represented  in  the  form: 

z(t,. )  =  h[x(t,. )  t,.  ]  +  v(t; )  (3-50) 

By  expanding  the  double-differenced  range  term  and  combing  terms  Equation  (3-49) 
becomes: 

(3-51) 

where  now 

=|(AVr^-hAV/;^+AV/«;^-hAV<)-t-AViV;^-AViV;^  (3-52) 

A 

The  term  modeled  as  a  white  noise,  is  the  combination  of  the  residual  error  terms  in 

the  double-differenced  equation  and  the  difference  in  the  ambiguity  terms.  As  explained 
previously,  even  though  the  individual  residual  error  terms 
(AVr“*,AV/“*,AVm''*,AVv^)  are  not  white  noises,  we  assume  the  combination  of  the 
four  terms  can  be  sufficiently  modeled  as  a  white  noise. 
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As  long  as  the  ambiguity  terms  in  Equation  (3-52)  and  (3-49)  cancel,  i.e. 

=  0,  modeling  as  a  white  noise  is  a  good  assumption.  However,  if 

the  two  ambiguity  terms  are  not  equal,  a  constant  bias  equal  to  the  difference  in  the 
ambiguity  terms  will  remain  in  the  ambiguity  corrected  double-differenced  carrier-phase 

measurements.  The  presence,  or  absence,  of  this  bias,  AVAT^  -  AVA^* ,  will  be  apparent 

in  the  elemental  filter’s  residuals.  This  bias  in  the  residuals  is  the  key  to  determining 
which  elemental  filter  contains  the  correct  ambiguity  set. 

Expanding  the  range  terms  and  substituting  state  variables  for  the  appropriate  terms 
as  described  in  Section  (3.3.2)  gives  the  following  equation: 

-^)^  +(/  -^2)'  +  (2“  -^3)']'' 

Again,  the  partial  derivative  matrix  H  must  be  calculated.  The  individual  partial 
derivatives  for  a  row  of  the  double-differenced  widelane  carrier-phase  measurements  are: 

mi)  ^  i(^''-^i)^+(3'“-^2)^+(2‘'-^3)^]^^ 

_j_| _ x”  -X, _ 1  (3-54) 

^[[(x'' -x.y+iy" -x,y  +{z‘’ -x,yr 


3l[x,tj] 

Sx^ 
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Sh[x,t]\ 
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These  can  be  combined  into  a  vector  representing  one  row  of  the  H  matrix  as: 


H“''=  -(e;-e;)  0  0  0  0  0  0 
A 


(3-57) 


where 


e“  =  [<?j  ^2  ^3]  =  unit  line  of  sight  vector  between  the  receiver  and  satellite  a. 

e*  =  [gj  ^2  63]  =  unit  line  of  sight  vector  between  the  receiver  and  satellite  b. 

The  entire  measurement  matrix  H  is  an  (n-l)-by-9  matrix  defined  by: 

-(e'-e^)  0  0  0  0  0  0 

A 

-(e'-e^)  0  0  0  0  0  0 

H  =  ^  (3-58) 

-(e’-e")  0  0  0  0  0  0 


The  measurement  error  covariance  matrix  for  the  elemental  filters  defined  by: 
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(3-59) 


Rfe)  for  ti=tj 
0  for  t-  ^  tj 


The  R  matrix  is  a  (n-l)-by-(n-l)  matrix  specified  by: 

R  = 

where 

=  double-differenced  carrier-phase  measurement  variance 
r  (,  =  double-differenced  carrier-phase  measurement  covariance 
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The  double-differenced  carrier-phase  measurement  variances  were  all  set  to  the  same 
value  as  were  the  double-differenced  carrier-phase  measurement  covariances.  The 
elemental  filter  measurement  error  values  used  in  this  research,  for  both  ground  and  flight 
test,  are  listed  in  Table  3-5.  Due  to  the  dynamics  associated  with  the  flight  test  data  the 
measurement  covariance  error  values  were  reduced  to  allow  the  conditional  probability  to 
more  rapidly  shift  between  elemental  filters. 


Table  3-5.  Elemental  Filter  Measurement  Covariance  Error  Values 


Term 

Definition 

Ground  Test  Value 

Flight  Test  Value 

^VA^'^VA^y 

Double-differenced  carrier- 
phase  variance  error 

(0.75m)^ 

(0.075m)^ 

Double-differenced  carrier- 
phase  covariance  error 

0.5*(0.75m)^ 

0.5*  (0.075m)^ 
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3.5.1.3  Discrete-Time  Models 

The  elemental  filter  model  equations  derived  above  must  be  converted  into  difference 
equations  in  order  to  be  implemented  easily  in  a  Kalman  filter.  Following  the  procedure 
described  in  Section  3.3.3,  again  assuming  =  the  elemental  filter  state 


transition  matrix  is  defined  as: 
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and  the  discrete  dynamics  driving  noise  matrix: 
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where 
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3.5.2  Conditional  Probability  Calculations 

The  conditional  probability  generation  block  shown  in  Figure  3-1  is  based  on 
Equations  (2-55)  and  (2-56),  derived  in  Section  2.4.2,  with  some  slight  modifications. 
Equations  (2-55)  and  (2-56)  are  repeated  here  for  clarity. 


where 


f s(<i  )|a.Z(rM )  h;  ’  ^i-1 )  Pj  (^1-1 ) 

A,(t,)  =  H,(r,)P,(tr)H/(f,)  +  R,(t,.) 

hih) = z,  -hk(trx^,] 


(2-55) 


(2-56) 
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The  success  of  the  MMAE  depends  on  the  conditional  probability  calculations  being 
able  to  distinguish  between  correct  and  incorrect  solutions,  i.e.  the  correct  ambiguity  set. 
The  conditional  probability  calculations  are  very  sensitive  to  differences  in  the  residuals 
of  the  elemental  filters.  Any  error,  or  noise,  in  the  state  estimates  and  covariances  of 
each  of  the  elemental  filters  increases  their  residual  values  and  degrades  the  overall 
performance  of  the  MMAE.  Therefore,  Equations  (2-55)  and  (2-56)  were  modified  to 
incorporate  the  latest  measurement;  i.e.  P^(rr)=>P^(r*)  and  Xi^(tT)=^  Using  the 

updated  state  estimate  and  covariance  values,  or  “post-fit”  residuals,  instead  of  the 
propagated  estimates  means  the  conditional  probability  density  function  in  Equation 
(2-56)  is  now  conditioned  on  the  time  history  of  measurements  up  to  and  including 
measurements  at  time  r,-.  The  post-fit  residuals  reduce  the  effect  of  propagation  errors 
within  the  elemental  filters  on  the  conditional  probability  calculations.  The  propagation 
errors  during  static  ground  tests  were  small  enough  that  the  algorithm  performed  well 
without  using  post-fit  residuals.  However,  during  flight  tests,  without  these  post-fit 
residuals,  the  propagation  errors  caused  the  measurement  residuals  to  be  large  enough 
that  the  conditional  probability  calculations  could  not  determine  which  elemental  filter 
contained  the  correct  ambiguity  set.  However,  careful  consideration  must  be  taken  before 
modifying  the  residual  calculations  so  as  not  to  disrupt  the  beneficial  qualities  of  an 
MMAE.  The  new  equations  using  the  post-fit  residuals  are  shown  below  in  Equations 
(3-63)  and  (3-64). 
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(3-63) 


where 


i=\  )|a,Z(j j )  (*i  I  )P  j  (^1-1 ) 


!((,  )|a,Z(l,- )  ) 


(2;r)'"'^|A,(r,.)| 


-1/2 


=exp{«} 


(3-64) 


A.(<, )  =  «.(», )P.(<,*)H/(r,)  +  R.((,) 
r,(l,)  =  z,-h[x, (I, •),<,] 


Recall  that  the  elemental  filters  are  identical,  except  for  the  candidate  ambiguity  set 
being  evaluated.  This  lead  to  the  value  of  |A^(r,.)|  being  the  same  for  all  elemental 
filters,  out  to  ten  significant  digits.  Therefore,  the  leading  coefficient, 

- — - -777 ,  was  removed  from  the  conditional  probability  calculations,  and  the 

conditional  probability  calculations  are  only  based  on  the  expj-^rj  (r,  )Ai(r,.)~‘r;t(t,  )|  in 

Equation  (2-56).  By  removing  the  leading  coefficient.  Equation  (2-56)  is  no  longer  a 
proper  conditional  probability  density,  and  tests  showed  there  was  no  significant  effect  on 
the  final  probability  weights,  other  than  simplifying  the  necessary  calculations. 


3.5.2.1  Lower  Probability  Bounds 

As  discussed  in  Section  2.4.2. 1,  a  conditional  probability  lower  bound  must  be  set  to 
ensure  the  MMAE  does  not  prematurely  ignore  information  from  a  particular  elemental 
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filter.  The  lower  probability  bound  (e)  for  each  elemental  filter  was  set  at  0.001.  A 
smaller  lower  probability  bound  (e  =  0.0005)  yielded  a  sluggish  response  in  the 
calculated  probabilities  due  to  changes  in  the  elemental  filters.  A  larger  lower  probability 
bound  (e  =  0.01)  reduced  the  accuracy  of  the  overall  MMAE  estimate  by  inappropriately 
forcing  large  weights  to  be  associated  with  incorrect  solutions.  The  lower  probability 
bound  of  E  =  0.001  proved  to  be  a  good  compromise  between  providing  a  quick  enough 
response  to  changes  in  the  MMAE  and  not  overly  weighting  the  MMAE  estimate  toward 
a  wrong  elemental  filter. 

3.5.3  Filter  Pruning 

An  important  aspect  of  the  overall  algorithm  not  depicted  in  Figure  3-1  is  the  concept 
of  filter  pruning.  Most  MMAEs  are  designed  to  look  constantly  for  a  change,  or  changes 
(failed  sensor  or  actuator,  etc.)  in  the  system  of  interest.  Because  the  MMAE  is 
constantly  looking  for  this  change,  each  elemental  filter  must  be  checked  at  every  sample 
time  to  determine  if  it  represents  the  correct  change  to  the  system.  However,  in  the 
integer  ambiguity  resolution  problem,  each  elemental  filter  is  used  to  check  a  constant 
bias,  or  ambiguity  set.  Because  integer  ambiguities  are,  by  definition,  constants  once  an 
ambiguity  set  is  determined  wrong,  it  will  always  be  wrong.  There  is  then  no  benefit  to 
keeping  that  particular  elemental  filter  within  the  MMAE. 

In  addition,  to  generate  an  accurate  navigation  solution,  the  MMAE  must  converge  to 
the  true  fixed-point  solution,  meaning  the  elemental  filter  based  on  the  correct  ambiguity 
set  must  absorb  the  entire  probability  weight.  However,  based  on  Equation  (2-57)  this 
cannot  happen  unless  there  is  only  one  elemental  filter  in  the  MMAE.  Therefore,  logic 
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was  added  to  the  algorithm  to  remove  an  elemental  filter  from  the  MMAE  if  the 
elemental  filter’s  probability  remained  at  or  below  the  lower  bound  e  for  more  than  a 
predetermined  number  of  consecutive  sample  periods.  The  number  of  consecutive 
sample  periods  was  determined  empirically  and  can  be  varied  for  specific  applications. 
For  ground  tests,  elemental  filters  were  removed  from  the  MMAE  after  10  consecutive 
sample  periods  below  the  lower  bound  (at  a  IHz  data  rate).  For  flight  tests,  elemental 
filters  were  removed  from  the  MMAE  after  40  consecutive  samples  below  the  lower 
bound  periods  (at  a  2Hz  data  rate).  The  additional  time  during  flight  test  was  necessary 
for  the  MMAE  to  determine  the  correct  ambiguity  set  due  to  the  larger  errors  caused  by 
the  dynamics  of  the  flight  tests.  The  overall  MMAE  state  estimate  is  then  simply  an 
application  of  Equation  (2-57)  on  the  elemental  filters  still  left  in  the  algorithm. 

By  adding  pruning  logic  to  the  overall  algorithm  design,  the  MMvAE  is  now  similar  to 
some  forms  of  moving  bank  MMAE  designs  [26,  27].  However,  instead  of  moving  a 
constant  number  of  elemental  filters  through  a  search  space  this  new  design  changes  the 
number  of  elemental  filters  used  to  examine  the  search  space. 

The  Bayesian  method  described  in  Equation  (2-57)  uses  all  of  the  remaining 
elemental  filters  to  provide  an  overall  state  estimate.  This  ensures  that  until  the  pruning 
logic  deletes  all  but  the  correct  elemental  filter,  the  filters  based  on  a  wrong  ambiguity  set 
will  have  some  affect  on  the  overall  MMAE  state  estimate.  Although  the  wrong 
elemental  filters  associated  conditional  probabilities  may  be  small,  their  solutions  can 
only  degrade  the  performance  of  the  MMAE. 

A  different  approach  would  be  to  use  the  maximum  a  posteriori,  or  MAP,  design  of 
the  MMAE.  The  MAP  design  uses  the  same  conditional  probability  calculations  as  a 
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normal  MMAE  design.  However,  in  the  MAP  design,  the  overall  MMAE  state  estimate 
is  the  state  estimate  of  the  elemental  filter  that  has  the  highest  conditional  probability.  If 
the  elemental  filter  with  the  highest  probability  is  the  correct  filter  then  the  MMAE 
estimate  is  not  degraded  by  the  solution  of  wrong  elemental  filters.  However,  this 
approach  ignores  the  information  available  in  the  other  elemental  filters.  In  addition,  if 
the  filter  with  the  highest  probability  is  the  wrong  filter,  as  it  occasionally  is  in  this 
research  as  shown  in  Section  4,  then  the  overall  state  estimate  may  be  worse  than  the 
Bayesian  approach. 

3.6  Summary 

This  chapter  discussed  the  overall  algorithm,  system  models,  and  techniques  used  in 
this  thesis.  The  floating-point  DGPS  Kalman  fitter  and  the  MMAE  structure  are  based  on 
previous  AFIT  research.  The  carrier-phase  ambiguity  set  generation  routine  was  a 
combination  of  the  LAMBDA  and  EASE  methods.  The  floating-point  DGPS  Kalman 
filter  model  and  the  elemental  filter  models  used  in  the  MMAE  were  presented  and 
described  in  detail.  In  addition,  changes  made  to  the  MMAE  theory  presented  in  Chapter 
2,  necessary  for  the  development  of  the  overall  algorithm,  were  described.  Chapter  4  will 
discuss  the  completed  results  based  on  both  ground  and  flight  tests. 
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4  Results  and  Analysis 


4.1  Overview 

This  chapter  presents  the  results  and  analysis  of  the  algorithm  developed  in  this 
thesis.  The  first  section  covers  four  static  ground  tests  cases,  and  the  second  section 
covers  the  three  flight  tests  cases.  The  results  presented  here  represent  the  best  possible 
tuning  within  the  scope  of  this  research.  Note  that,  the  purpose  of  this  research  is  to 
analyze  the  use  of  a  MMAE  structure  to  determine  carrier-phase  ambiguities  and  generate 
an  accurate  GPS  navigation  solution,  not  to  design  a  production  quality  flight-worthy 
system. 

4.2  Ground  Tests 

Ground  tests  were  conducted  at  Wright-Patterson  Air  Force  Base  (AFB)  Ohio,  from 
1  October  99  to  31  October  99.  Results  are  presented  for  four  different  cases,  with 
baseline  distances  ranging  from  250  m  to  100  km.  Dual  frequency  Ashtech  Z-Surveyor 
receivers  were  used  to  collect  data  at  both  the  reference  and  remote  sites.  The  reference 
station  was  located  at  a  surveyed  site  on  the  roof  of  the  Air  Force  Institute  of  Technology. 
The  remote  site  was  a  tripod-mounted  antenna  and  receiver  placed  at  various  distances 
from  the  reference  station. 
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4.2.1  Case  Definitions 

Each  ground  test  case  consisted  of  approximately  30  minutes  of  data  collected  at  a  1- 
Hz  rate.  The  data  included  LI  C/A  pseudorange  and  carrier-phase  measurements,  and  L2 
pseudorange  (semi-codeless)  and  carrier-phase  measurements.  In  all  cases  the  lower 
probability  bound  for  each  elemental  filter  was  0.001,  as  justified  in  Section  3. 5.2.1,  and 
each  had  an  equal  initial  probability  weighting.  The  pruning  logic  was  set  to  delete 
elemental  filters  if  their  conditional  probability  remained  below  the  lower  limit  for  more 
than  10  consecutive  sample  periods.  The  floating-point  and  elemental  filters  were  tuned 
in  accordance  with  the  values  specified  in  Chapter  3.  The  MMAE  used  between  15  and 
100  filters,  and  was  manually  picked  depending  on  the  case.  More  elemental  filters  were 
used  as  the  expected  difficultly  of  resolving  the  carrier-phase  ambiguities  increased,  due 
to  increased  baseline  distances.  Eventually,  the  algorithm  would  automatically  pick  the 
number  of  elemental  filters  based  on  an  expected  level  of  difficulty. 

The  first  two  ground  test  cases  were  chosen  to  demonstrate  the  ability  of  the  MMAE 
to  “instantaneously”  resolve  the  carrier-phase  ambiguities,  using  either  widelane  or  LI 
carrier-phase  measurements.  The  last  two  cases  were  chosen  to  demonstrate  the  ability  of 
the  MMAE  to  correctly  resolve  the  carrier-phase  ambiguities  at  medium  and  long 
baseline  distances.  Ground  test  case  specifics  are  shown  below  in  Table  4-1. 


Table  4-1.  Ground  Test  Case  Parameter  Description 


Case 

Baseline 

Distance 

Carrier-Phase 
Measurement  Type 

Elemental 

Filters 

Data 

Rate 

Length  of 
Data  Set 

1 

250  m 

Widelane 

15 

IHz 

~  30  Minutes 

2 

LI 

25 

~  30  Minutes 

3 

Widelane 

100 

■EH 

~  30  Minutes 

4 

■UliIrlM 

Widelane 

100 

■EH 

~  30  Minutes 
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For  the  ground  tests,  the  truth  position  and  the  true  carrier-phase  ambiguities  were 
determined  by  the  Ashtech  Office  Suite  software  [1].  The  position  errors  are  defined  as: 


^error  ^filter  ^true 


(4-1) 


where  the  subscript  “filter”  refers  to  either  the  floating-point  filter  or  the  MMAE 
solution,  as  appropriate,  and  the  subscript  “true”  refers  to  the  Ashtech  Office  Suite 

solution.  Note  that  is  not  the  absolute  truth  but  it  is  the  solution  the  MMAE  is 
striving  to  achieve,  namely  the  fixed-integer  solution. 


4.2.2  Case  1: 250  m  Baseline,  Widelane  Measurements 

Case  1  is  designed  to  demonstrate  the  near  instantaneous  ambiguity  resolution 
capability  of  the  MMAE.  The  time  history  position  performance  of  the  floating-point 
filter  and  the  MMAE,  for  the  entire  data  run,  are  shown  in  Figure  4-1  and  Figure  4-2, 
respectively.  In  both  figures,  the  solid  line  represents  the  computed  position  error  as 
described  by  Equation  (4-1).  The  dashed  line  represents  the  filter-predicted  standard 
deviation,  not  the  actual  standard  deviation  of  the  results.  At  the  beginning  of  the  data 
run,  the  MMAE  was  initialized  to  the  floating-point  filter,  which  had  initial  position 
errors  between  0.2-0.7  meters.  As  shown  in  Figure  4-2,  the  MMAE  corrected  these 
errors  to  less  than  0.1  meters  in  one  sample  period.  The  steps  in  the  filter-predicted 
standard  deviation  within  the  first  100  seconds  are  due  to  the  algorithm  pruning  bad 
elemental  filters  as  described  in  Section  3.5.3.  The  filter-predicted  standard  deviations 
for  the  Up  direction  in  Figure  4-2  are  at  approximately  ±0.6m,  just  outside  the  plot  limits. 
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Pos  error  (m) 


The  time  averaged  mean  position  error  and  standard  deviation  of  the  floating-point 
filter  solution  and  the  MMAE  solution  are  listed  below  in  Table  4-2.  The  results  are 
represented  in  the  geodetic  (East,  North,  Up  (ENU))  reference  frame.  These  results  show 
that  the  MMAE  solution  was  an  order  of  magnitude  more  accurate  than  the  floating-point 
filter  solution. 


Table  4-2.  Case  1  -  Position  Solution  Error  Comparison 


East 

Non 

th 

_ Ufi _ 1 

Mean 
Error  (m) 

Std  (m) 

Mean 
Error  (m) 

Std  (m) 

Mean 
Error  (m) 

Std  (m) 

Floating-Point 

Filter 

-0.107 

0.118 

-0.059 

0.138 

-0.172 

0.205 

MMAE 

0.023 

0.019 

-0.024 

0.010 

0.059 

0.042 

The  first  minute  of  the  conditional  probability  time  histories  for  two  of  the  15 
elemental  filters  are  shown  in  Figure  4-3.  In  Figure  4-3,  elemental  filter  4  (EF  4) 
contains  the  correct  carrier-phase  ambiguity  set,  and  elemental  filter  3  (EF  3)  contains  an 
incorrect  carrier-phase  ambiguity  set.  After  the  first  probability  calculation,  elemental 
filter  4  absorbs  approximately  98%  of  the  conditional  probability,  essentially 
instantaneously  resolving  the  ambiguity.  The  remaining  2%  of  the  conditional 
probability  is  distributed  among  the  other  14  incorrect  elemental  filters.  The  conditional 
probability  of  12  of  the  14  incorrect  elemental  filters  is  small  enough  that  the  filters 
conditional  probability  had  to  be  artificially  set  at  the  lower  bound  of  0.001.  Ten  seconds 
later  the  filter  pruning  logic  deletes  12  of  the  14  incorrect  elemental  filters  (those  whose 
probability  stayed  at  the  lower  bound)  and  elemental  filter  4  absorbs  approximately 
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EF  4  EF  3 


99.8%  of  the  conditional  probability.  Within  25  seconds,  the  filter  pruning  logic  deletes 
the  remaining  2  incorrect  elemental  filters.  Elemental  filter  4  is  the  only  remaining 
elemental  filter  in  the  MMAE,  and  because  it  contains  the  correct  ambiguity  set,  the 
MMAE  converges  to  the  true  fixed-integer  solution. 

1 

0.5 

0 
1 
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GPS  Time  (sec) 

Local  Time  (hh:mm:ss) 

Figure  4-3.  Case  1  -  Sample  Elemental  Filter  Conditional  Probability  Time 

Histories,  First  Minute 

4.2.3  Case  2: 250  m  Baseline,  LI  Measurements 

Case  2  is  identical  to  Case  1,  except  the  algorithm  uses  LI  carrier-phase 
measurements  instead  of  widelane  measurements.  The  wavelength  of  an  LI  cycle  (~19 
cm)  is  approximately  one-fourth  the  wavelength  of  a  widelane  measurement  (~86  cm). 
This  increases  the  number  of  candidate  ambiguity  sets  within  a  given  search  space, 
making  the  correct  ambiguity  set  harder  to  resolve.  However,  once  resolved,  an  LI 
fixed-integer  solution  is  more  accurate  than  a  widelane  fixed-integer  solution,  due  to 
amplification  of  multipath  and  ionospheric  errors  in  the  widelane  measurements. 
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The  performance  of  the  floating-point  filter  and  the  MMAE  are  shown  in  Figure  4-4 
and  Figure  4-5,  respectively.  Again,  the  solid  line  represents  the  computed  position  error 
and  the  dashed  line  represents  the  filter-predicted  standard  deviation. 

In  this  case,  the  floating-point  filter  initializes  with  position  errors  of  less  than  0.2 
meters,  and  the  MMAE  corrects  these  errors  to  less  than  0.03  meters  within  one  sample 
period.  Note  that  the  vertical  scale  in  Figure  4-5  is  five  times  smaller  than  in  Figure  4-4. 

The  floating-point  filter  uses  the  same  tuning  values  for  this  case  as  Case  1,  even 
thought  different  carrier-phase  measurements  are  used.  This  might  explain  why  the 
floating-point  filter  solutions  are  consistently  outside  of  the  filter  predicted  covariance 
bounds.  For  this  case  the  floating-point  filter  performance  is  satisfactory,  however,  an 
operational  algorithm  would  use  separate  tuning  values  for  different  measurement  types. 


GPS  Time  (sec) 
Local  Time  (hh:mm) 


Figure  4-4.  Case  2  -  Floating-Point  Filter  Position  Errors 
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Figure  4-5.  Case  2  -  MMAE  Position  Errors 


The  time  averaged  mean  position  error  and  standard  deviation  of  the  floating-point 


filter  solution  and  the  MMAE  solution  (ENU  frame)  are  shown  in  Table  4-3.  As 


expected,  the  results  indicate  that  the  MMAE  solution  is  at  least  an  order  of  magnitude 


more  accurate  than  the  floating-point  filter  solution,  and  the  LI  solution  is  more  accurate 


than  the  wi delane  solution  shown  in  Case  1. 


Table  4-3.  Case  2  -  Position  Solution  Error  Comparison 
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EF  24  EF  21 

Rcsidusl  3  Rssidusl  1  Cond  Prob  R©sidu3l  3  Rosidu©!  1  Cond  Prob 


The  conditional  probability  time  history,  and  two  of  the  measurement  residual  (post- 
fit  residual)  time  histories  for  elemental  filters  21  and  24  are  shown  below  in  Figure  4-6 
and  Figure  4-7,  respectively.  Only  the  first  two  minutes  of  the  data  run  are  shown. 
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Figure  4-6.  Case  2  -  Elemental  Filter  21  Probability  and  Residual  Time  History 

(Correct  Ambiguities) 
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Figure  4-7.  Case  2  -  Elemental  Filter  24  Probability  and  Residual  Time  History 

(Incorrect  Ambiguities) 


4-9 


Both  of  the  residuals  shown  in  Figure  4-6  for  elemental  filter  21  appear  almost  zero- 
mean  and  white;  characteristics  of  a  residual  in  a  properly  designed  and  tuned  Kalman 
filter.  However,  in  Figure  4-7  the  first  residual  for  elemental  filter  24  appears  to  be 
growing  more  negative,  while  measurement  residual  3  shows  a  small  but  noticeable  bias. 
These  small  differences  in  residuals  are  enough  to  cause  the  conditional  probability  to 
shift  rapidly  toward  elemental  filter  number  21.  The  straight  horizontal  line  in  the 
elemental  filter  21  residual  plots,  approximately  25  seconds  after  the  start  of  the  data  run, 
are  due  to  the  pruning  logic,  which  deleted  the  elemental  filter  from  the  MMAE  at  that 
point. 

Because  of  the  small  value  of  its  residuals,  elemental  filter  21  (the  correct  elemental 
filter)  absorbs  approximately  97%  of  the  probability  during  the  first  conditional 
probability  calculation.  The  remaining  3%  is  divided  among  the  other  24  incorrect 
elemental  filters.  In  a  manner  very  similar  to  Case  1,  the  pruning  logic  starts  deleting 
elemental  filters  with  probabilities  at  the  lower  bound  for  10  consecutive  samples. 
Within  40  seconds,  the  pruning  logic  deletes  all  of  the  incorrect  elemental  filters,  and  the 
MMAE  converges  to  the  true  fixed-integer  solution. 

The  results  from  the  first  two  cases  show  that  for  short  baselines,  an  MMAE  almost 
instantly  determines  which  elemental  filter  contains  the  correct  carrier-phase  ambiguity 
set.  In  addition,  the  MMAE  solution  is  less  erratic  and  at  least  and  order  of  magnitude 
more  accurate  than  the  DGPS  floating-point  filter  solution. 
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4.2.4  Case  3: 28.1  km  Baseline,  Widelane  Measurements 

Case  3  is  designed  to  demonstrate  the  performance  of  the  MMAE  using  a  medium 
baseline  distance.  The  performance  of  the  floating-point  filter  and  the  MMAE  are  shown 
in  Figure  4-8  and  Figure  4-9.  The  floating-point  filter  (and  hence  the  MMAE)  initializes 
with  a  much  larger  than  expected  position  error  (~  0.5-2.5  meters),  due  to  residual 
differential  errors  not  canceled  in  the  floating-point  filter.  These  residual  differential 
errors  are  a  problem,  because  they  are  used  by  the  FASF  routine  to  generate  candidate 
ambiguity  sets,  and  they  decrease  the  ability  of  the  MMAE  to  determine  which  elemental 
filter  contains  the  correct  ambiguity  set.  Consequently,  the  correct  elemental  filter  is  not 
instantly  obvious  to  the  MMAE  in  the  conditional  probability  calculations.  For  the  first  3 
minutes  of  the  data  run,  the  conditional  probability  shifts  between  an  incorrect  elemental 
filter  and  the  correct  elemental  filter,  causing  the  transients  in  the  MMAE  position  errors 
seen  in  Figure  4-9.  Even  thought  initially  the  MMAE  heavily  weights  an  incorrect 
elemental  filter,  its  position  solution  accuracy  is  comparable  to  the  floating-point  filter 
solution  during  the  initial  transients. 

The  other  advantage  of  the  MMAE  solution  is  that  it  transitions  very  smoothly  from  a 
floating-point  filter  like  solution  to  the  true  fixed-integer  solution.  There  are  no 
discontinuities,  or  jumps,  in  its  position  solution  (typical  in  other  ambiguity  resolution 
routines)  as  it  approaches  the  fixed-integer  solution.  The  MMAE  is  reset  after  it  lost  lock 
on  a  satellite  (at  ~19;01),  and  even  thought  the  transients  are  not  as  large  as  the  initial 
transients,  the  MMAE  again  transitions  very  rapidly  and  smoothly  from  the  floating-point 
solution  to  the  fixed-integer  solution. 
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Figure  4-8.  Case  3  -  Floating-Point  Position  Errors 
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Figure  4-9.  Case  3  -  MMAE  Position  Errors 
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As  in  the  previous  cases,  position  error  comparisons  between  floating-point  filter 
solution  and  the  MMAE  solution  are  shown  in  Table  4-4.  Again,  as  expected,  the 
MMAE  solution  is  significantly  more  accurate  than  the  floating-point  solution,  even  with 
the  large  transient  errors  averaged  into  the  overall  MMAE  accuracy.  The  large  initial 
transient  errors  are  also  the  cause  of  the  larger  than  expected  standard  deviation  in  the 
MMAE  solution  listed  in  Table  4-4. 


Table  4-4.  Case  3  -  Position  Solution  Error  Comparison 


East 

Non 

th 

_ Up _ 1 

Mean 
Error  (m) 

Std  (m) 

Mean 
Error  (m) 

Std  (m) 

Mean 
Error  (m) 

Std  (m) 

Floating-Point 

Filter 

0.355 

0.074 

0.112 

0.159 

0.533 

0.327 

MMAE 

0.018 

0.151 

0.060 

0.098 

0.205 

0.500 

The  conditional  probability  time  histories  for  the  four  elemental  filters  of  interest  in 
this  case  are  shown  Figure  4-10.  Initially  elemental  filter  68  (EF  68)  contains  the  correct 
carrier-phase  ambiguity  set,  and  elemental  filter  5  (EF  5)  contains  an  incorrect  carrier- 
phase  ambiguity  set.  At  approximately  19:01  local  time,  the  mobile  receiver  lost  lock  on 
a  low  elevation  satellite.  At  this  time  the  algorithm  reset  all  of  the  elemental  filters  in  the 
MMAE  to  the  floating  filter  solution,  generated  new  candidate  ambiguity  sets  for  each 
elemental  filter,  and  reset  the  elemental  filter  probabilities  to  equal  values.  After  being 
reset  elemental  filter  77  (EF  77)  contains  the  correct  ambiguity  set,  and  elemental  filter 
58  (EF  58)  contains  an  incorrect  ambiguity  set.  Resetting  the  probabilities  explains  the 
instant  probability  drop  in  elemental  filter  68,  and  the  probability  spikes  in  elemental 
filters  58  and  77. 
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Note  that  resetting  the  MMAE  was  done  for  progranaming  simplicity,  and  is  not 
characteristic  of  an  operational  algorithm.  Completely  reinitializing  the  MMAE  removes 
any  information  or  advantage  already  gained  by  the  MMAE.  The  performance  of  the 
algorithm  could  be  improved  by  using  the  information  already  available  in  the  MMAE  to 
reset  the  elemental  filter  and  more  efficiently  handle  changes  in  the  satellite  constellation. 
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Figure  4-10.  Case  3  -  Sample  Elemental  Filter  Conditional  Probability  Time 

Histories 


The  first  six  minutes  of  conditional  probability  time  history,  and  four  representative 
measurement  residual  time  histories  for  elemental  filters  5  and  68  are  shown  in  Figure 
4-11  and  Figure  4-12,  respectively.  As  shown  in  Figure  4-10  and  Figure  4-11,  the 
MMAE  initially  shifts  the  probability  weight  to  elemental  filter  5.  Even  though  two  of 
the  four  measurement  residuals  (residuals  2  and  3)  are  relatively  large,  the  other  two 
residuals  are  small  enough  to  cause  this  initial  probability  shift  toward  elemental  filter  5. 
However,  the  2"^  and  3’’''  residuals  in  elemental  filter  5  are  slowly  increasing  in 
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magnitude.  This  increase  causes  the  conditional  probability  to  shift  away  from  elemental 
filter  5.  Eventually,  the  elemental  filter  5  residuals  are  consistently  large,  which  reduces 
the  probability  to  near  zero  and  causes  the  pruning  logic  to  remove  the  filter  from  the 
MMAE.  Because  the  elemental  filter  68  residuals  are  much  closer  to  zero  mean  and 
white,  the  MMAE  quickly  starts  to  weight  elemental  filter  68  (which  contains  the  correct 
ambiguity  set).  In  approximately  3.5  minutes,  the  pruning  logic  deletes  the  incorrect 
elemental  filters  and  the  MMAE  converges  to  the  correct  filter  and  the  true  fixed-integer 
solution. 
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Figure  4-11.  Case  3  -  Elemental  Filter  5  Probability  and  Residual  Time  History 
(Incorrect  Ambiguities),  First  Six  Minutes 


4-15 


345000  345060  345120  345180  345240  345300  345360 

18:50  18:51  18:52  18:53  18:54  18:55  18:56 


GPS  Time  (sec) 

Local  Time  (hh:mm) 

Figure  4-12.  Case  3  -  Elemental  Filter  68  Probability  and  Residual  Time  History 
(Correct  Ambiguities),  First  Six  Minutes 
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4.2.5  Case  4: 100  km  Baseline,  Widelane  Measurements 

The  purpose  of  Case  4  is  to  demonstrate  the  ability  of  the  MMAE  to  determine  the 
correct  ambiguity  set  under  the  more  difficult  conditions  associated  with  long  a  baseline 
distance.  The  baseline  distance  was  increased  to  -100  km  to  increase  the  affect  of  the 
residual  differential  errors  in  the  ambiguity  set  generation  routine  and  the  conditional 
probability  calculations. 

The  initialization  errors  from  the  floating-point  filter  for  this  case  were  large  like 
Case  3,  and  had  a  major  effect  on  the  candidate  ambiguity  sets  generated  by  the  FASF 
subroutine.  In  fact,  none  of  the  100  candidate  ambiguity  sets  initially  generated  by  the 
FASF  routine  were  the  correct  ambiguity  set.  This  meant  the  MMAE  could  not  converge 
to  the  correct  solution,  as  the  correct  solution  was  not  in  one  of  the  elemental  filters. 

Similar  to  Case  3,  the  transients  in  the  MMAE  position  solution  in  Figure  4-14  are 
due  to  the  initial  high  weighting  of  an  incorrect  elemental  filter  (elemental  filter  3  as 
shown  in  Figure  4-15).  However,  elemental  filter  3  does  not  contain  the  correct 
ambiguity  set,  and  the  MMAE  quickly  shifts  the  probability  weight  to  elemental  filter  62. 
The  MMAE  converges  to  elemental  filter  62  in  approximately  6  minutes,  even  though 
elemental  filter  62  does  not  contain  the  correct  ambiguity  set.  Recall  that  none  of  the 
elemental  filters  contain  the  correct  ambiguity  set. 

The  MMAE  is  reset  at  approximately  22:06:30  local  time  when  the  remote  receiver 
losses  lock  on  a  low  elevation  satellite.  At  this  time,  the  correct  ambiguity  set  is  in  one  of 
the  100  elemental  filters  (filter  67)  and  the  MMAE  converges  to  this  filter  and  the  correct 
solution  within  90  seconds  of  being  reset. 
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Figure  4-15.  Case  4  -  Sample  Elemental  Filter  Conditional  Probability  Time 

Histories 

The  true  ambiguities  and  the  ambiguities  used  by  the  elemental  filters  are  shown 
below  in  Table  4-5.  Elemental  filter  3  contains  three  incorrect  ambiguities  (each  one 
cycle  off  the  correct  value),  while  elemental  filter  62  only  contains  one  incorrect 
ambiguity.  Although  it  contains  an  incorrect  ambiguity  set,  elemental  filter  62  was  the 
best  ambiguity  set  available,  among  the  100  elemental  filters  (remember,  the  correct 
ambiguity  set  was  not  included  in  the  100  elemental  filters).  The  ambiguity  set  in 
elemental  filter  62  is  considered  the  best  set,  as  it  is  the  ambiguity  set  closest  to  the  true 
ambiguity  set.  The  other  elemental  filters  contain  ambiguity  sets  with  more  than  one 
incorrect  ambiguity  or  an  ambiguity  set  with  one  incorrect  ambiguity  more  than  one  cycle 
off  the  correct  value.  The  difference  between  the  ambiguity  set  in  an  elemental  filter  and 
the  true  ambiguity  set  determines  the  size  of  the  measurement  bias  in  the  ambiguity 
corrected  double-differenced  carrier-phase  measurements  described  in  Section  3. 5. 1.2.  In 
general,  the  larger  the  difference  between  the  two  ambiguity  sets,  the  larger  the  error  in 
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the  computed  position  solution.  The  underlined  and  italicized  entries  are  the  incorrect 
ambiguity  values,  and  the  N/A  for  elemental  filter  67  is  due  to  a  satellite  being  lost  from 
view,  and  therefore  that  measurement  (and  associated  ambiguity  term)  being  lost. 


Table  4-5.  Case  4  -  True  Ambiguities  vs.  Elemental  Filter  Ambiguities 


True 

Ambiguity 

1094483 

2268285 

3812363 

4213446 

4210673 

3709768 

7338024 

EF#3 

1094483 

226cS286 

3812363 

4213446 

4210672 

3709768 

7338023 

EF#62 

1094483 

2268285 

3812363 

4213446 

4210673 

3709768 

7338023 

EF#67 

1094483 

2268285 

3812363 

4213446 

4210673 

3709768 

N/A 

The  time  averaged  mean  error  and  standard  deviations  of  the  floating  filter  compared 
to  those  of  the  MMAE  design  are  shown  in  Table  4-6.  Again,  the  MMAE  performance 
produced  much  better  results  than  the  floating-point  filter,  even  though  the  MMAE  could 
not  initially  converge  to  the  true  fixed-integer  solution. 


Table  4-6.  Case  4  -  Position  Solution  Error  Comparison 


East 

North 

_ Up _ 1 

Mean 
Error  (m) 

Std  (m) 

Mean 
Error  (m) 

Std  (m) 

Mean 
Error  (m) 

Std  (m) 

Floating-Point 

Filter 

0.167 

0.223 

-0.545 

0.400 

0.437 

0.192 

MMAE 

0.014 

0.207 

-0.055 

0.167 

0.253 

0.121 

This  case  clearly  demonstrates  some  of  the  strengths  of  using  an  MMAE  to  resolve 
carrier-phase  ambiguities  and  in  the  process  generate  a  navigation  position  solution. 
Figure  4-14  clearly  shows  the  MMAE  solution  transitioning  very  smoothly  from  a 
floating-point  filter  like  solution  to  a  fixed-integer  solution.  Even  though  the  true  fixed- 
integer  solution  is  not  reached  (it  is  not  initially  available)  the  MMAE  transitions  to  the 


4-20 


best  solution  available.  Also  notice,  that  the  MMAE  fixed-integer  solution  is  more 
accurate  and  more  stable  than  the  floating-point  solution  even  though  it  is  using  and 
incorrect  ambiguity  set.  After  the  MMAE  is  reset,  due  to  the  loss  of  a  satellite  from  view, 
the  correct  ambiguity  set  is  available  and  the  MMAE  again  transitions  very  quickly  and 
smoothly  to  the  true  fixed-integer  solution. 

None  of  the  initial  candidate  ambiguity  sets  generated  by  the  algorithm  contain  the 
correct  ambiguity  set.  As  described  in  Section  3.4  the  initial  candidate  ambiguity  sets  are 
generated  after  only  one  cycle  of  the  floating-point  filter.  After  one  cycle  the  floating¬ 
point  filter  may  not  have  an  accurate  estimate  of  the  carrier-phase  ambiguities.  However, 
after  the  MMAE  is  reset  the  correct  ambiguity  set  was  generated  by  the  Z-transform  and 
EASE  subroutines.  By  this  time  the  floating-point  filter  has  a  more  accurate  estimate  of 
the  carrier-phase  ambiguities,  which  increases  the  probability  of  the  correct  ambiguity  set 
being  included  in  one  of  the  elemental  filters.  Allowing  the  floating-point  filter  a  few 
propagate  and  update  cycles  before  using  its  ambiguity  estimates  to  generate  the 
candidate  ambiguity  sets  may  improve  the  performance  of  the  carrier-phase  ambiguity  set 
generation  routine  and  the  overall  algorithm. 
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4.3  Flight  Tests 

Flight  tests  were  conducted  at  Edwards  AFB  CA,  from  5  October  2000  to  17  October 
2000.  Appendix  A  describes  the  complete  series  of  flight  tests.  Out  of  the  twelve  flight 
conditions  outlined  in  Appendix  A,  results  are  presented  for  three  different  cases.  Due  to 
flight  limitations,  data  was  only  collected  for  the  primary  data  points  described  in 
Appendix  A. 

The  reference  receiver  station  was  a  dual  frequency  Ashtech  Z-Surveyor  receiver 
located  at  a  surveyed  site  at  the  Air  Force  Flight  Test  Center  (AFFTC)  Range  Division, 
Edwards  AFB,  CA.  Data  presented  in  this  section  is  represented  in  a  geodetic  (East, 
North,  Up  (ENU))  reference  frame  with  the  origin  located  at  the  reference  receiver 
location.  The  coordinates  of  the  reference  receiver  station  are  listed  in  Appendix  A.  The 
test  aircraft  was  a  T-38A  aircraft,  modified  to  carrier  an  Ashtech  Z-surveyor  receiver  in 
the  rear  cockpit. 

4.3.1  Case  Definitions 

Each  flight  test  case  consisted  of  4  to  20  minutes  of  data  collected  at  a  2  Hz  rate  as 
outlined  in  Appendix  A.  For  each  case,  the  test  aircraft  was  flown  in  an  arc  at  a  constant 
speed,  distance  and  altitude  from  the  reference  receiver.  LI  C/A  -  code  pseudorange  and 
carrier-phase  measurements  and  L2  pseudorange  (semi-codeless)  and  carrier-phase 
measurements  were  collected  by  both  the  reference  receiver  and  the  receiver  onboard  the 
test  aircraft.  The  reference  receiver  was  at  an  altitude  of  -690  m  while  the  aircraft  was 
flown  at  an  altitude  of  at  least  3,400  m.  In  all  cases,  the  lower  probability  bound  for  each 
elemental  filter  was  set  at  0.001  and  each  had  an  equal  initial  probability  weighting.  The 
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pruning  logic  in  the  algorithm  was  set  to  delete  elemental  filters  if  their  conditional 
probability  remained  below  the  lower  limit  for  more  than  40  consecutive  sample  periods 
(20  seconds).  The  pruning  time  was  increased  to  allow  the  MMAE  more  time  to 
determine  which  elemental  filter  contained  the  correct  ambiguity  set,  and  was  necessary 
due  to  the  additional  errors  added  by  the  aircraft  dynamics,  not  present  during  the  ground 
tests.  The  floating-point  filter  and  the  elemental  filters  were  tuned  in  accordance  with  the 
values  specified  in  Chapter  3.  For  all  flight  test  cases,  the  MMAE  within  the  algorithm 
used  50  elemental  filters.  The  number  of  elemental  filters  was  chosen  based  on  the 
expected  difficultly  of  resolving  the  carrier-phase  ambiguities.  However,  as  discussed  in 
the  next  section,  the  carrier-phase  ambiguity  set  generation  subroutine  did  not  perform 
well  under  flight  test  conditions. 

The  first  flight  test  case  was  chosen  to  demonstrate  the  performance  of  the  MMAE 
using  a  nominal  flight  test  profile.  The  second  case  uses  the  same  flight  profile  but 
demonstrates  the  problems  associated  with  removing  elemental  filters  from  an  MMAE. 
The  last  case  was  chosen  to  demonstrate  the  ability  of  the  MMAE  to  correctly  resolve  the 
carrier-phase  ambiguities  at  an  increased  baseline  distance.  The  flight  test  case  specifics 
are  shown  below  in  Table  4-7. 


Table  4-7.  Flight  Test  Case  Parameter  Description 


Case 

Baseline 

Distance 

Pressure 

Altitude 

Carrier-Phase 

Measurement 

Type 

Elemental 

Filters 

Data 

Rate 

Length  of 
Data  Set 

1 

16.1  km 
(10  naut.  mi) 

3,400  m 
(10,000  ft) 

Widelane 

50 

2Hz 

~  4  minutes 

2 

16.1  km 
(10  naut.  mi) 

3,400  m 
(10,000  ft) 

Widelane 

50 

2Hz 

~  4  minutes 

3 

32.2  km 
(20  naut.  mi) 

3,400  m 
(10,000  ft) 

Widelane 

50 

2Hz 

~  6  minutes 
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For  the  flight  tests,  the  true  carrier-phase  ambiguities  were  determined  by  the  Ashtech 
Office  Suite  software  [1].  The  truth  position  was  then  calculated  using  a  fixed-integer 
carrier-phase  DGPS  routine,  using  the  true  ambiguities  determined  by  Ashtech  Office 
Suite  software.  For  all  flight  test  cases,  the  position  errors  are  defined  as: 

^ error  ~  ^filter  ~  ^true  (4-2) 

where  the  subscript  “filter”  refers  to  either  the  floating-point  filter  or  the  MMAE 
solution,  and  the  subscript  “true”  refers  to  the  fixed-integer  calculated  truth  position. 
The  “truth”  position  generated  by  using  the  correct  carrier-phase  ambiguities  and  a  fixed- 
integer  carrier-phase  DGPS  routine  is  not  an  absolute  truth  position  solution,  but  it  is  the 
solution  the  MMAE  is  striving  to  achieve. 

The  position  errors  described  above  are  only  used  to  generate  the  position  error  plots 
in  the  following  sections.  They  are  not  used  to  calculate  time  averaged  mean  errors  and 
standard  deviations.  Because  of  the  way  the  “truth”  is  generated  for  the  flight  tests,  once 
the  MMAE  converges  to  the  fixed-integer  solution  the  two  solutions  (“truth”  and 
MMAE),  are  virtually  identical,  and  it  is  impossible  to  determine  if  the  errors  are  due  to 
the  MMAE  solution  or  the  truth  solution.  In  addition,  the  small  amount  of  data  used  in 
each  flight  test  case,  along  with  the  initial  transients  of  the  MMAE  solution,  limits  the 
usefulness  of  a  time  averaged  mean  error.  Therefore,  no  time  averaged  mean  error  or 
standard  deviations  are  presented  for  the  flight  test  cases. 
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4.3. 1.1  Ambiguity  Set  Generation 

For  the  flight  test  cases  analyzed,  none  of  the  50  candidate  ambiguity  sets  initially 
generated  by  the  carrier-phase  ambiguity  set  generation  routine  (Z-transform,  FASF  and 
inverse  Z-transforms  sub-routines)  were  the  correct  ambiguity  set.  This  meant  the 
MMAE  could  not  converge  to  the  correct  solution,  as  the  correct  solution  was  not  in  one 
of  the  elemental  filters.  Although  not  changed  for  this  research,  allowing  the  floating¬ 
point  filter  a  few  propagate  and  update  cycles  before  using  its  ambiguity  estimates  to 
generate  the  candidate  ambiguity  sets  may  improve  the  performance  of  the  carrier-phase 
ambiguity  set  generation  routine. 

The  main  emphasis  of  this  research  is  to  investigate  the  ability  of  an  MMAE  to 
distinguish  the  correct  ambiguity  set;  therefore,  the  correct  ambiguity  set  was  artificially 
placed  in  the  first  elemental  filter.  In  addition,  an  ambiguity  set  with  one  incorrect 
ambiguity  (one  cycle  off  the  correct  ambiguity)  was  artificially  placed  in  the  second 
elemental  filter. 

Placing  these  ambiguity  sets  in  the  elemental  filters  allowed  the  performance  of  the 
MMAE  to  be  analyzed  in  two  ways.  First,  the  ability  of  the  MMAE  to  determine  which 
ambiguity  set  is  correct  could  be  analyzed.  Second,  the  ability  of  the  MMAE  to 
distinguish  between  closely  related  solutions  could  also  be  determined. 

The  above  modification  was  implemented  to  validate  the  “proof  of  concept” 
algorithm  developed  in  this  thesis.  An  operational  system  would  require  advanced 
routines  to  increase  the  probability  that  the  correct  ambiguity  set  is  in  one  of  the 
elemental  filters. 
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4.3.2  Case  5:  16.1  km  Baseline,  10,000  ft  Altitude,  Widelane  Measurements 

The  purpose  of  this  flight  test  case  is  to  demonstrate  the  performance  of  the  MMAE 
using  a  nominal  flight  test  profile.  The  time  history  position  performance  of  the  floating¬ 
point  filter  and  the  MMAE,  for  the  entire  data  run,  are  shown  in  Figure  4-16  and  Figure 
4-17  respectively.  The  errors  are  represented  in  a  geodetic  ENU  reference  frame  with  the 
origin  at  the  reference  receiver  location.  In  both  figures,  the  solid  line  represents  the 
computed  position  error  as  described  in  Equation  (4-2).  The  dashed  line  represents  the 
filter-predicted  standard  deviation,  and  not  the  actual  standard  deviation  of  the  errors. 

As  discussed  in  Section  4.3. 1.1,  the  correct  ambiguity  set  was  not  generated  by  the 
FASF  subroutine,  so  the  correct  ambiguity  set  was  artificially  placed  in  elemental  filter  1. 
The  conditional  probability  time  histories  for  elemental  filters  1,  2, 14,  and  41  are  shown 
in  Figure  4-18.  The  conditional  probability  shifts  back  and  forth  between  elemental 
filters  14  and  41  for  the  first  45  seconds  of  the  data  run.  At  approximately  45  seconds, 
elemental  filter  1  rapidly  absorbs  the  conditional  probability  and  within  75  seconds,  the 
pruning  logic  deleted  all  of  the  other  elemental  filters.  The  probability  shifts  between 
elemental  filters  is  the  reason  for  the  large  errors  initially  present  in  the  MMAE  position 
solution.  However,  once  the  MMAE  converges  to  the  correct  elemental  filter,  the 
position  solution  is  very  stable  and  within  a  few  centimeters  of  the  true  position  solution. 

The  large  “spike”  in  the  floating-point  filter  solution  at  8:07:00  is  due  to  a  bad 
pseudorange  measurement.  Because  the  MMAE  only  uses  carrier  phase  measurements, 
the  bad  pseudorange  measurement  only  affected  the  floating-point  filter  solution. 
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Figure  4-18.  Case  5  -  Sample  Elemental  Filter  Conditional  Probability  Time 

Histories 

The  trae  ambiguities  and  the  ambiguities  used  by  the  elemental  filters  are  shown  in 
Table  4-8.  The  underlined  and  italicized  values  represent  incorrect  ambiguity  values. 
Elemental  filter  1  contains  the  correct  ambiguity  set,  while  elemental  filter  2  contains 
only  one  incorrect  ambiguity.  Elemental  filter  14,  which  twice  manages  to  absorb  over 
95%  of  the  conditional  probably,  contains  only  one  correct  ambiguity.  Elemental  filter 
14  even  contains  one  ambiguity  17  cycles  off  the  correct  value.  Finally,  elemental  filter 
41  does  not  contain  any  correct  ambiguities.  The  ambiguity  sets  contained  in  the 
elemental  filter  are  much  further  off  of  the  true  ambiguity  set  than  expected.  This  causes 
the  initial  transients  to  have  much  larger  errors  than  in  the  ground  test  cases.  The 
ambiguity  sets  are  a  direct  result  of  allowing  the  floating-point  filter  only  one  propagate 
and  update  cycle  before  its  estimates  are  used  to  generate  candidate  ambiguity  sets. 
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Table  4-8.  Case  5  -  True  Ambiguities  vs.  Elemental  Filter  Ambiguities 


True 

Ambiguity 

-2320781 

-1205221 

210188 

1885413 

12121638 

5587938 

EF#1 

-2320781 

-1205221 

210188 

1885413 

12121638 

5587938 

EF#2 

-2320781 

-1205221 

210188 

1885413 

12121639 

5587938 

EF#14 

-2320798 

-1205218 

210188 

1885417 

12121629 

5587916 

EF#41 

-2320792 

-1205219 

210191 

1885425 

12121636 

5587925 

The  conditional  probability  time  histories  and  three  of  the  measurement  residual  time 
histories  for  elemental  filters  1  and  14  are  shown  in  Figure  4-19  and  Figure  4-20, 
respectively. 
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Local  Time  (hh:mm:ss) 

Figure  4-19.  Case  5  -  Elemental  Filter  1  Probability  and  Residual  Time  History 

(Correct  Ambiguities) 
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Figure  4-20.  Case  5  -  Elemental  Filter  14  Probability  and  Residual  Time  History 

(Incorrect  Ambiguities) 

For  the  first  10  to  15  seconds  of  the  data  run,  none  of  the  elemental  filter  residuals  are 
small  enough  for  the  conditional  probability  calculations  to  determine  which  elemental 
filter  contains  the  correct  ambiguity  set.  Then  the  value  of  the  2““*  and  4*  residuals  in 
elemental  filter  14  become  very  small,  even  though  it  contains  an  ambiguity  set  with 
ambiguities  far  from  the  correct  values.  In  fact,  the  elemental  filter  14  residuals  are 
smaller  than  the  same  residuals  in  elemental  filter  1,  demonstrating  how  residual  DGPS 
errors  can  make  it  difficult  to  distinguish  the  true  ambiguity  set.  After  approximately  45 
seconds,  the  first  residual  in  elemental  filter  14  starts  to  increase  in  value.  It  is  this 
increase  (along  with  an  increase  in  the  residuals  not  shown)  that  causes  the  conditional 
probability  weight  to  shift  away  from  elemental  filter  14  and  toward  elemental  filter  1. 
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Although  elemental  filter  2  contained  only  one  incorrect  ambiguity,  its  residuals  (not 
shown)  were  initially  larger  than  the  residuals  of  elemental  filter  1, 14  or  41.  The  pruning 
logic  quickly  deleted  elemental  filter  2  from  the  MMAE,  because  its  probability  was 
below  the  lower  bound  for  the  first  20  seconds. 

Although  the  MMAE  transitions  smoothly  to  the  true  fixed-integer  solution  it 
emphasizes  the  affect  of  a  poorly  performing  candidate  ambiguity  set  generation  routine 
on  the  overall  performance  of  the  algorithm.  The  candidate  ambiguity  generation  routine 
did  not  generate  ambiguity  sets  close  to  the  true  ambiguity  set,  and  for  demonstration 
purposes  the  correct  set  was  artificially  placed  in  the  MMAE.  Because  the  incorrect 
candidate  ambiguity  sets  are  so  far  from  the  correct  values,  they  cause  much  larger  errors 
and  more  rapid  position  shifts  in  the  overall  solution  than  would  incorrect  ambiguity  sets 
closer  to  the  true  ambiguity  set.  This  case  also  emphasizes  the  problem  of  rapidly 
shifting  the  probability  to  an  elemental  filter  whose  probabilities  look  “good”  for  a  short 
period  of  time.  Tuning  the  MMAE  to  slow  down  the  conditional  probability  shifts,  or 
adding  additional  tests  on  the  elemental  filter  residuals,  such  as  a  whiteness  test  [15,  25], 
may  also  improve  the  overall  performance  of  the  algorithm. 
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4.3.3  Case  6:  16.1  km  Baseline,  10,000 ft  Altitude,  Widelane  Measurements 

Case  6  is  a  separate  data  run,  but  it  is  the  same  flight  profile  as  Case  5.  This  case, 
although  not  specifically  designed  too,  demonstrates  the  problems  associated  with 
removing  elemental  filters  from  and  MMAE.  The  time  history  position  performance  of 
the  floating-point  filter  and  the  MMAE,  for  the  entire  data  run,  are  shown  in  Figure  4-21 
and  Figure  4-22,  respectively.  In  both  figures,  the  solid  line  represents  the  computed 
position  error  and  the  dashed  line  represents  the  filter  predicted  standard  deviation. 
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Figure  4-21.  Case  6  -  Floating-Point  Filter  Position  Errors 
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Figure  4*22.  Case  6  -  MMAE  Position  Errors 

As  seen  in  the  previous  cases,  the  initial  transients  in  the  solution  are  due  to  the 
MMAE  shifting  the  probability  weight  to  an  incorrect  elemental  filter.  However,  unlike 
previous  cases,  the  MMAE  converges  to  a  solution  with  a  significant  bias  from  the  true 
solution.  This  large  bias  in  the  position  solution  of  the  MMAE  is  because  the  MMAE 
converged  to  the  wrong  elemental  filter,  even  thought  the  correct  ambiguity  set  is 
available  in  elemental  filter  1. 

The  conditional  probability  time  histories  for  elemental  filters  1,  2,  34,  38  and  45  are 
shown  in  Figure  4-23.  The  conditional  probability  calculations  never  shift  the  probability 
weight  toward  elemental  filter  1.  Elemental  filter  1  is  deleted  from  the  MMAE  23 
seconds  from  the  start  of  the  data  run.  Because  the  correct  ambiguity  set  is  no  longer 
available,  the  MMAE  can  never  converge  to  the  correct  fixed-integer  solution.  The 
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conditional  probability  shifts  between  elemental  filters  45  and  34  before  being 
completely  shifted  to  elemental  filter  38.  The  MMAE  converges  to  elemental  filter  38, 
and  the  wrong  solution,  53  seconds  after  the  start  of  the  data  run. 
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Figure  4-23.  Case  6  -  Conditional  Probability  Time  Histories 


The  true  ambiguities  and  the  ambiguities  used  by  the  elemental  filters  are  shown  in 
Table  4-9.  Elemental  filter  1  contains  the  correct  ambiguity  set,  while  elemental  filter  38, 
the  filter  the  MMAE  converges  to,  has  four  incorrect  ambiguities,  one  ten  cycles  off  the 
correct  value. 
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Table  4-9.  Case  6  -  True  Ambiguities  vs.  Elemental  Filter  Ambiguities 


True 

Ambiguity 

-2320781 

-1205221 

210188 

-3181105 

12121638 

5587938 

EF#1 

-2320781 

-1205221 

210188 

-3181105 

12121638 

5587938 

EF#2 

-2320781 

-1205221 

210188 

-3181105 

12121639 

5587938 

EF#34 

-2320780 

-1205228 

210193 

-3181092 

12121635 

5587947 

EF#38 

-2320786 

-1205221 

210191 

-3181095 

12121638 

5587932 

EF#45 

-2320793 

-1205234 

210191 

-3181109 

12121614 

5587943 

The  probability  time  history  and  three  of  the  measurement  residual  time  histories  for 
elemental  filters  1,  38  and  45  are  shown  below  in  Figure  4-24,  Figure  4-25,  and  Figure 
4-26  respectively.  As  shown  in  Figure  4-26,  the  elemental  filter  45  measurement 
residuals  at  the  start  of  the  data  run,  are  small  enough  that  the  MMAE  shifts  the 
conditional  probability  to  elemental  filter  45.  However,  this  shift  is  only  for  a  few 
seconds,  because  the  residuals  rapidly  increase  in  value,  causing  the  MMAE  to  shift  the 
probability  to  elemental  filter  34. 

The  elemental  filter  1  measurement  residuals  are  large  enough  at  the  beginning  of  the 
data  run  to  shift  the  conditional  probability  away  from  the  filter.  Measurement  residual  3 
in  Figure  4-24  is  actually  getting  larger  at  the  start  of  the  data  run.  In  contrast,  two  of  the 
elemental  filter  38  residuals  (residuals  2  and  4)  are  very  close  to  zero.  Approximately  20 
seconds  after  the  start  of  the  data  run  all  six  (even  the  three  not  shown)  of  the  elemental 
filter  38  residuals  are  very  close  to  zero.  This  causes  the  rapid  shift  in  probability  toward 
the  filter.  None  of  the  other  elemental  filter’s  residuals  are  small  enough  to  cause  the 
MMAE  to  shift  the  probability  away  from  elemental  filter  38.  The  pruning  logic 
eventually  deletes  the  other  elemental  filters  and  the  MMAE  is  forced  to  converge  to 
elemental  filter  38. 
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Figure  4-24.  Case  6  -  Elemental  Filter  1  Probability  and  Residual  Time  Histories 

(Correct  Ambiguities) 
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Figure  4-25.  Case  6  -  Elemental  Filter  38  Probability  and  Residual  Time  Histories 

(Incorrect  Ambiguities) 
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Figure  4-26.  Case  6  -  Elemental  Filter  45  Probability  and  Residual  Time  Histories 

(Incorrect  Ambiguities) 

This  is  the  first  case  where  the  MMAE  converges  to  the  wrong  elemental  filter  (and 
ambiguity  set),  even  though  the  correct  ambiguity  set  is  in  one  of  the  elemental  filters. 
Because  an  incorrect  filter’s  residuals  are  smaller  than  the  correct  filter’s  residuals  for  a 
period  of  time,  the  MMAE  shifts  the  probability  toward  the  wrong  filter.  The  pruning 
logic  then  deletes  the  correct  filter  from  the  solution  before  the  MMAE  has  a  chance  to 
shift  the  probability  back  to  the  correct  elemental  filter.  Close  to  the  end  of  the  data  run 
all  of  the  elemental  filter  38  measurement  residuals  are  increasing.  This  would  normally 
cause  the  MMAE  to  shift  the  probability  toward  another  elemental  filter,  possibly  the 
correct  elemental  filter.  However,  because  the  pruning  logic  removed  all  of  the  other 
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elemental  filters,  the  MMAE  must  continue  to  use  elemental  filter  38  as  the  correct 
solution. 

This  case  clearly  shows  the  weakness  in  the  current  algorithm’s  pruning  logic  design. 
Once  the  pruning  logic  deletes  the  elemental  filters  from  the  MMAE  there  is  no  way  to 
add  additional  elemental  filters  to  check  or  recheck  different  ambiguity  sets.  An  MMAE 
is  a  Multiple  Model  Adaptive  Estimator.  Removing  elemental  filters,  in  a  desire  to 
achieve  the  true  fixed-integer  solution,  limits  the  ability  of  the  MMAE  to  adapt  or  change 
its  solution.  When  the  MMAE  is  reduced  to  one  elemental  filter  it  can  no  longer  adapt  to 
any  changes,  as  it  is  forced  to  use  the  one  remaining  hypothesis,  even  if  it  is  the  incorrect 
hypothesis.  These  problems  point  to  many  possible  improvements  in  the  current 
algorithm.  A  simple  modification  of  not  deleting  the  elemental  filters  from  the  MMAE 
would  allow  the  MMAE  to  continually  adapt  to  changes.  The  overall  MMAE  estimate 
can  then  be  a  blending  of  only  those  elemental  filters  with  a  high  conditional  probability 
of  being  the  correct  solution.  Additional  modifications  might  include  deleting  “bad” 
elemental  filters,  but  then  adding  elemental  filters  to  recheck  ambiguity  sets,  or  ambiguity 
sets  close  to  the  elemental  filter  that  has  the  highest  probability.  Additional 
improvements  to  limit  the  chances  of  the  MMAE  converging  the  wrong  solution  might  be 
to  continuously  cross  check  the  MMAE  solution  with  the  floating-point  solution.  The 
two  solutions  should  be  closely  related  (within  a  few  meters),  and  if  they  are  not  the 
MMAE  could  be  reinitialized,  or  different  ambiguity  sets  could  be  used. 
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4.3.4  Case  7:  32.2  km  Baseline,  10,000 ft  Altitude,  Widelane  Measurements 

Case  7  is  similar  to  Cases  5  and  6,  except  the  baseline  distance  to  the  test  aircraft  was 
increased  to  32.2  km.  The  time  history  position  performance  of  the  floating-point  filter 
and  the  MMAE  for  the  entire  data  run  are  shown  in  Figure  4-27  and  Figure  4-28, 
respectively.  As  shown  in  Figure  4-28,  the  MMAE  very  quickly  and  smoothly  converges 

to  the  correct  elemental  filter  and  the  true  fixed-integer  solution. 

The  large  “spike”  in  the  floating-point  filter  solution  at  7:47:30  is  due  to  a  bad 
pseudorange  measurement.  Again,  because  the  MMAE  only  uses  carrier  phase 
measurements,  the  bad  pseudorange  measurement  only  affected  the  floating-point  filter 

solution. 

Once  again  the  correct  ambiguity  set  was  placed  in  elemental  filter  1,  and  an  incorrect 
ambiguity  set  (one  incorrect  ambiguity)  was  placed  in  elemental  filter  2.  The  true 
ambiguities  and  the  ambiguities  in  the  elemental  filters  are  shown  below  in  Table  4-10. 


Table  4-10.  Case  7  -  True  Ambiguities  vs.  Elemental  Filter  Ambiguities 


True 

Ambisuity 

1081928 

-1205221 

210188 

-1885413 

5587938 

1081928 

-120522in 

210188  n 

-3181105^ 

5587938 

1081928 

-1205221 

210188 

-3181105 

5587939 

4-39 


The  conditional  probability  time  histories  and  three  of  the  measurement  residual  time 
histories  for  elemental  filters  1  and  2  are  shown  in  Figure  4-29  and  Figure  4-30 
respectively.  Figure  4-30  shows  the  elemental  filter  2  residuals  for  the  entire  data  run, 
even  though  the  elemental  filter  was  deleted  from  the  MMAE  21  seconds  after  the  start  of 
the  data  run. 

As  seen  in  Figure  4-29,  the  elemental  filter  1  residuals  were  all  very  close  to  zero  for 
the  entire  data  run.  This  allowed  elemental  filter  1  to  absorb  over  90  percent  of  the 
conditional  probability  in  the  first  10  seconds.  The  MMAE  converged  to  elemental  filter 
1  (the  correct  solution)  32  seconds  after  the  start  of  the  data  run.  In  contrast,  the 
elemental  filter  2  residuals  (Figure  4-30)  either  show  a  constant  bias  or  and  increase  in 
magnitude  over  the  7  minutes  of  data. 

The  MMAE  converged  very  quickly  to  the  correct  solution  due  to  the  almost  ideal 
behavior  of  the  elemental  filter  1  residuals.  The  performance  of  this  algorithm  (and  most 
other  carrier-phase  ambiguity  resolution  algorithms)  is  very  dependent  on  the  “quality”  of 
the  data  set.  A  clean  data  set,  or  even  a  period  of  clean  data  within  a  data  set,  can 
dramatically  improve  the  performance  of  an  algorithm,  and  is  a  very  typical  problem  with 
most  carrier-phase  ambiguity  resolution  routines. 
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Figure  4-29.  Case  7  -  Elemental  Filter  1  Probability  and  Residual  Time  Histories 

(Correct  Ambiguities) 
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Figure  4-30.  Case  7  -  Elemental  Filter  2  Probability  and  Residual  Time  Histories 

(Incorrect  Ambiguities) 
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4.4  Summary 

This  chapter  presented  the  performance  results  of  an  MMAE-based  carrier-phase 
integer  ambiguity  resolution  algorithm.  The  results  were  broken  up  into  ground  and 
flight  test  cases,  with  each  case  representing  a  difference  distance,  or  distance  and 
altitude  from  the  reference  receiver.  Each  case  used  real  world  data  from  Ashtech  Z- 
surveyor  receivers  at  both  the  reference  and  remote  (or  airborne)  locations. 

During  ground  tests,  the  MMAE  was  very  efficient  at  determining  the  carrier-phase 
integer  ambiguities  at  short  baseline  distances.  As  expected,  the  MMAE  required  more 
time  to  determine  the  correct  integer  ambiguity  set  as  the  baseline  distances  increased.  In 
addition,  even  though  the  accuracy  of  the  new  algorithm  degraded  as  baseline  distances 
increased,  the  new  algorithm  provided  a  more  accurate  position  solution  than  a  floating¬ 
point  DGPS  solution.  In  most  cases,  the  MMAE  solution  was  as  accurate,  if  not  more 
accurate,  than  the  floating-point  solution  as  it  converged  to  the  true  fixed-integer  solution. 
In  all  cases,  the  MMAE  solution  was  accurate  to  the  centimeter  level  once  it  converged  to 
the  correct  elemental  filter.  During  the  final  ground  test,  the  new  algorithm  demonstrated 
the  ability  of  the  MMAE  to  converge  to  the  best  possible  ambiguity  set  in  the  elemental 
filters  even  if  the  correct  ambiguity  set  is  not  available.  Even  thought  the  MMAE 
converged  to  the  best  possible  solution,  improvements  to  algorithm  to  increase  the 
probability  of  the  correct  solution  being  available  within  the  MMAE  were  suggested. 

During  flight  tests,  the  correct  integer  ambiguity  set  was  not  initially  generated  by  the 
candidate  ambiguity  set  generation  subroutine,  which  suggests  numerous  improvements 
are  required  in  the  candidate  ambiguity  set  generation  subroutine.  Therefore,  the  correct 
carrier-phase  integer  ambiguity  set  was  artificially  placed  in  one  of  the  elemental  filters. 
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Once  the  correct  solution  was  available  to  the  MMAE,  the  MMAE  was  generally  able  to 
determine  the  correct  carrier-phase  integer  ambiguities,  but  the  results  were  inconsistent. 
In  two  of  the  three  flight  test  cases  presented,  the  MMAE  converged  to  the  correct 
solution.  In  the  remaining  case,  the  MMAE  converged  to  the  wrong  solution  even 
thought  the  correct  solution  was  available.  In  this  case,  the  correct  filter  was  prematurely 
removed  from  the  MMAE  by  the  pruning  logic  before  the  MMAE  was  able  to  converged 
to  it.  This  problem  suggests  improvements  not  in  the  MMAE,  but  in  the  pruning  logic. 
By  removing  the  elemental  filters  the  MMAE  cannot  adapt  to  changes,  and  is  forced  to 
converge  to  whatever  solution  is  available.  Numerous  simple  changes  may  be  made  in 
the  pruning  logic  to  improve  the  overall  performance  of  the  algorithm. 

The  actual  position  solution  accuracy  of  the  new  algorithm,  under  flight  test 
conditions,  was  not  always  more  accurate  than  the  floating-point  filter.  Before  the 
MMAE  converged  to  the  correct  elemental  filter,  when  the  conditional  probability  shifted 
to  an  incorrect  elemental  filter,  the  MMAE  solution  was  less  accurate  than  the  floating¬ 
point  filter  solution.  These  large  initial  errors  were  again  due  to  the  poor  candidate 
ambiguity  sets  generated  initially  used  in  the  elemental  filters,  again  suggesting 
improvements  in  the  floating-point  filter  and  its  role  in  the  candidate  ambiguity  set 
generation  subroutine.  However,  once  the  MMAE  converged  to  the  correct  solution,  it 
generated  a  navigation  solution  with  centimeter  level  accuracy 

In  all  ground  test  and  flight  test  cases,  even  when  the  MMAE  converged  to  the  wrong 
solution,  it  very  smoothly  and  quickly  converged  to  a  fixed-integer  solution.  There  were 
no  discontinuities  or  jumps  in  its  solution  typical  of  other  carrier-phase  ambiguity 
resolution  routines. 
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5  Conclusions  and  Recommendations 


5.1  Overview 

This  research  presented  the  theory,  modeling  methodology,  and  results  of  a  GPS 
carrier-phase  integer  ambiguity  resolution  algorithm  based  on  a  multiple  model  Kalman 
filter.  Previous  research  in  this  area  considered  the  use  of  MMAE  algorithms  for  many 
purposes  such  as  detection  of  failed  actuators  or  sensors  and  GPS  spoofing  and  jamming. 
Other  research  included  extensive  study  into  stand-alone  integer  ambiguity  resolution 
algorithms,  including  both  ambiguity  set  generation  and  ambiguity  set  determination. 
This  research  represents  the  first  known  use  of  an  MMAE  algorithm  in  the  carrier-phase 
ambiguity  resolution  problem.  The  focus  was  not  in  the  generation  of  ambiguity  sets,  but 
in  the  determination  of  the  correct  ambiguity  set  while  using  all  available  information  to 
improve  the  GPS  navigation  solution. 

Two  different  extended  Kalman  filter  models  were  developed  for  this  thesis:  a 
floating-point  DGPS  filter  and  the  elemental  filters  within  the  MMAE  structure.  Both  of 
the  filters  incorporated  a  linear  dynamics  model  with  a  nonlinear  measurement  model. 
The  floating-point  DGPS  filter  was  used  to  initialize  the  MMAE,  but  more  importantly, 
its  outputs  were  used  to  generate  the  candidate  ambiguity  set  used  in  the  elemental  filters 
in  the  MMAE.  The  elemental  filters  in  the  MMAE  were  each  based  on  a  different 
candidate  ambiguity  set  in  order  to  determine  the  correct  integer  ambiguity  set,  while 
using  all  available  information  to  provide  an  accurate  navigation  solution.  The  standard 
MMAE  structure  was  also  modified  to  incorporate  filter  pruning  and  to  use  the  latest 
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measurement,  or  “post-fit”  residuals,  in  the  MMAE  conditional  probability  calculations, 
rather  than  conventional  filter  residuals.  Both  static  ground  tests  and  flight  tests  were 
conducted  to  verify  the  performance  of  the  new  algorithm. 

5.2  Conclusions 

The  new  algorithm  developed  in  this  thesis  performed  extremely  well  during  the 
ground  tests  and  flight  tests  analyzed  throughout  this  research.  The  new  algorithm  was 
able  to  determine  the  correct  carrier-phase  ambiguities  quickly  and  provided  an  accurate 
navigation  solution.  The  new  algorithm  showed  the  additional  benefit  of  transitioning 
rapidly  and  smoothly  from  a  floating-point  position  solution  to  a  fixed-integer  position 
solution.  In  all  cases,  once  the  MMAE  converged  to  the  correct  elemental  filter,  the 
solution  was  an  order  of  magnitude  more  accurate  than  the  floating-point  filter  solution, 
and  provided  a  centimeter  level  accurate  solution. 

Of  the  ground  test  cases  presented  in  Chapter  4,  each  shows  the  new  algorithm  will 
converge  to  the  correct  integer  ambiguity  set,  and  hence  the  correct  fixed-integer 
solution,  as  long  as  the  correct  ambiguity  set  is  in  one  of  the  elemental  filters.  The  new 
algorithm  essentially  instantaneously  resolved  the  correct  carrier-phase  ambiguities  for 
short  baseline  distances.  As  expected,  the  algorithm’s  performance  degraded  as  the 
baseline  distance  from  the  reference  receiver  was  increased.  The  MMAE  also  showed 
the  ability  to  converge  to  the  best  possible  solution,  if  the  correct  solution  was  not 
available  as  a  hypothesis  within  the  MMAE.  In  all  cases,  the  MMAE  smoothly 
approached  the  fixed-integer  solution  without  the  discontinuities  typical  of  most  carrier- 
phase  ambiguity  resolution  algorithms.  In  addition,  while  transitioning  to  the  fixed- 
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integer  solution,  the  MMAE  solution  was  generally  more  accurate  than  the  floating-point 
solution. 

For  flight  tests,  the  algorithm  showed  performance  improvements  over  current 
algorithms  even  though  the  results  were  not  as  consistent  as  during  the  ground  tests.  For 
flight  tests,  the  new  algorithm  was  again  able  to  determine  the  correct  carrier-phase 
integer  ambiguities  (after  the  correct  solution  was  artificially  placed  in  an  elemental 
filter),  but  the  results  were  inconsistent  due  to  problems  in  the  pruning  logic.  In  two  of 
the  three  cases  presented,  the  MMAE  converged  to  the  correct  solution.  Before  the 
MMAE  converged  to  the  correct  elemental  filter,  when  the  conditional  probability  shifted 
to  an  incorrect  elemental  filter,  the  MMAE  solution  was  sometimes  less  accurate  than  the 
floating-point  filter  solution.  However,  these  results  were  a  direct  function  of  the 
ambiguity  sets  generated  by  the  floating-point  filter  and  the  ambiguity  set  generation 
subroutine,  and  not  a  function  of  the  performance  of  the  MMAE.  Improving  the 
performance  of  the  floating-point  filter  and  the  ambiguity  set  generation  subroutine  will 
improve  the  performance  of  the  overall  algorithm.  The  MMAE  solution  smoothly  and 
rapidly  transitioned  to  a  fixed-integer  solution,  and  once  the  MMAE  converged  to  the 
correct  solution,  its  accuracy  was  at  the  cm-level  —  an  order  of  magnitude  more  accurate 
than  the  floating-point  solution. 

In  the  one  remaining  flight  test  case,  the  MMAE  converged  to  the  wrong  solution, 
even  thought  the  correct  solution  was  available  in  an  elemental  filter.  The  MMAE 
initially  shifted  the  probability  weight  to  an  incorrect  filter,  and  then  the  pruning  logic 
deleted  the  correct  filter  before  the  conditional  probability  could  be  shifted  back  to  the 
correct  elemental  filter.  The  MMAE  converged  to  the  wrong  solution  due  to  faulty 
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pruning  logic  and  not  as  a  result  of  the  MMAE  design.  This  case  illustrated  the  need  for 
careful  engineering  judgment  in  removing  elemental  filters  from  the  MMAE.  All  of  the 
tradeoffs  must  be  considered  before  removing  elemental  filters  and  reducing  or 
eliminating  the  ability  of  an  MMAE  to  adapt  to  changes  in  the  overall  system.  Numerous 
simple  changes  in  the  eurrent  pruning  logic  will  greatly  improve  the  performance  of  the 
algorithm. 

Although  not  the  focus  of  this  research,  the  ability  of  the  new  algorithm  to  resolve  the 
carrier-phase  ambiguities  correctly  depends  directly  on  the  correct  ambiguity  set  being 
hj^othesized  by  one  of  the  elemental  filters.  The  current  ambiguity  set  generation 
algorithm  does  not  guarantee  that  the  correct  ambiguity  set  is  available  in  one  of  the 
elemental  filters.  While  this  thesis  focused  on  the  ability  of  an  MMAE  to  determine  the 
correct  carrier-phase  ambiguity  set,  any  future  implementation  of  the  algorithm  created  in 
this  research  must  include  advanced  routines  to  improve  the  probability  of  the  correct 
ambiguity  being  available  as  a  solution  in  the  MMAE. 

The  ground  test  and  flight  test  profiles  used  in  this  research  were  adequate  to  test  the 
performance  of  the  algorithm  and  determine  the  effects  of  error  sources  in  a  DGPS 
application.  This  research  has  shown  that  under  ground  and  flight  test  conditions,  an 
MMAE  can  quickly  and  accurately  resolve  GPS  carrier-phase  ambiguities.  For  all  but 
one  case,  the  algorithm  converged  to  the  correct  ambiguity  set  and  provided  a  centimeter 
level  accurate  navigation  solution.  The  limited  results  presented  in  this  thesis  are  highly 
encouraging  and  warrant  further  research. 
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5.3  Recommendations 


This  initial  attempt  at  using  a  multiple  model  adaptive  filter  to  determine  the  correct 
carrier-phase  integer  ambiguities  performed  well  enough  to  continue  development. 
However,  there  are  many  modifications  that  may  be  made  to  the  current  algorithm  to 
improve  its  performance. 

The  following  recommendations  are  provided  to  extend  this  research  and  increase  the 
usefulness  of  the  current  algorithm: 

1.  The  performance  of  the  MMAE  directly  depends  on  the  correct  ambiguity  set  being 
hypothesized  by  one  of  the  elemental  filters.  To  improve  the  performance  of  the 
algorithm  it  should  be  modified  to: 

□  Incorporate  a  larger  number  of  elemental  filters  to  improve  the  chance  of  the 
correct  ambiguity  set  being  hypothesized  by  one  of  the  elemental  filters. 
However,  if  the  current  Bayesian  approach  (weighted  blending  of  the 
elemental  filter  solutions)  is  used  to  generate  the  overall  MMAE  state 
estimate,  increasing  the  number  of  elemental  filters  will  increase  the  effect  of 
incorrect  solutions  on  the  final  MMAE  estimate.  In  addition,  adding 
elemental  filters  increases  the  required  computational  power  and  memory. 

□  Allow  the  floating-point  filter  more  time  to  estimate  the  carrier-phase 
ambiguities  before  its  estimate  is  used  to  generate  candidate  ambiguity  sets. 
Currently,  the  candidate  ambiguity  sets  are  generated  after  only  one 
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propagate  and  update  cycle  of  the  floating-point  filter.  Allowing  the 
floating-point  filter  more  than  one  cycle  before  using  its  solution  to  generate 
the  candidate  ambiguity  sets  would  improve  the  likelihood  of  the  correct 
ambiguity  set  being  generated  by  the  Z-transform  and  FASF  routines. 

□  Reset  some  elemental  filters  with  ambiguity  sets  “close  to”  (as  determined  by 
the  appropriate  search  space)  the  “correct”  elemental  filter  once  the  MMAE 
has  converged  to  an  elemental  filter.  If  the  “correct”  elemental  filter  contains 
the  true  integer  ambiguity  set,  the  MMAE  will  again  converge  to  it  as  the 
fixed-integer  solution. 

□  Reset  the  elemental  filters  every  M  sample  periods  using  the  current 
candidate  ambiguity  set  generation  subroutine.  This  method  is  similar  to  the 
previous  recommendation,  but  alleviates  the  need  to  determine  what 
ambiguity  sets  are  “close  to”  the  “correct”  elemental  filter.  It  also  allows  the 
candidate  ambiguity  set  generation  subroutine  to  use  the  most  recent 
ambiguity  estimates  from  the  floating-point  filter. 

□  Perform  some  type  of  ‘likelihood  test”  [25]  on  the  elemental  filters  to 
provide  insight  into  the  “correctness”  of  the  solution  and  the  ambiguity  sets 
in  the  elemental  filters.  If  the  solution  has  a  low  likelihood  of  being  correct, 
the  elemental  filter(s)  may  be  reset  with  different  ambiguity  sets. 
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2.  The  current  algorithm  design  uses  a  Bayesian  estimator  (weighted  blending  of  the 
elemental  filter  solutions)  for  the  final  MMAE  state  estimates.  To  achieve  the  true 
fixed-integer  solution  with  the  Bayesian  estimator,  the  algorithm  also  uses  pruning 
logic  to  delete  “bad”  elemental  filters.  To  improve  the  performance  of  the  algorithm, 
the  final  MMAE  estimate,  and  avoid  the  problem  of  possibly  pruning  the  correct 
elemental  filter,  the  algorithm  should  be  modified  to: 

□  Not  prune  bad  elemental  filters,  or  not  prune  all  but  one  elemental  filter  from  the 
MMAE.  This  recommendation  removes  the  problem  of  the  MMAE  losing  its 
ability  to  adapt  to  system  changes,  and  the  possibility  of  removing  the  correct 
elemental  filter  from  the  MMAE.  The  overall  MMAE  solution  would  be  a  type  of 
blending  or  weighted  blending  of  the  elemental  filters  remaining  in  the  MMAE. 

□  Change  the  Bayesian  estimator  to  a  blended  weighting  scheme  similar  to  the  one 
described  in  reference  [10],  where  the  overall  estimate  is  still  a  weighted  average 
of  elemental  filter  solutions.  However,  only  those  elemental  filters  with  a 
conditional  probability  higher  than  some  predetermined  threshold  are  used  in  the 
final  MMAE  estimate.  If  the  threshold  is  set  high  enough  (but  not  restrictively 
high)  this  recommendation  removes  the  need  to  eliminate  elemental  filters,  and  if 
the  conditional  probability  of  the  “bad”  elemental  filters  is  below  the  threshold, 
then  the  final  MMAE  estimate  will  be  the  true  fixed-integer  solution. 
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3.  The  floating-point  filter  position  solution  is  only  used  to  initialize  and  reinitialize  the 
MMAE.  The  two  solutions  should  be  closely  related,  generally  within  a  few  meters 
based  on  the  accuracy  of  the  floating-point  filter.  Therefore,  the  floating-point  filter 
could  be  used  to  cross  check  the  MMAE  solution.  If  the  solutions  differ  by  a 
predetermined  amount,  the  MMAE  could  be  reinitialized,  or  different  ambiguity  sets 
could  be  generated  and  used  in  the  MMAE 

4.  Modify  the  MMAE  design  to  also  include  information  from  a  “whiteness”  test  [15] 
[25]  when  trying  to  determine  which  elemental  filter  contains  the  correct  carrier- 
phase  ambiguity  set.  A  whiteness  test  may  not  be  an  efficient  way  to  initially 
determine  the  correct  carrier-phase  ambiguity,  but  after  a  few  sample  periods  it  may 
add  enough  additional  information  to  stop  the  MMAE  from  incorrectly  weighting  an 
elemental  filter  whose  residuals  are  momentarily  small  as  they  pass  with  a  slope 
through  zero. 

5.  Increase  the  order  of  the  floating-point  filter  to  include  states  for  residual  differential 
tropospheric  errors  and  multipath  errors  to  improve  the  performance  of  the  floating¬ 
point  filter.  Improving  the  performance  of  the  floating-point  filter  will  improve  the 
performance  of  the  candidate  ambiguity  subroutine,  and  hence  increase  the  likelihood 
that  an  elemental  filter  will  contain  the  correct  ambiguity  set.  In  addition,  a  better 
floating-point  filter  solution  will  help  the  initialization  and  therefore  the  transient 
response  of  the  MMAE. 
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6.  Incorporate  logic  to  handle  the  addition  of  new  satellites  to  the  GPS  solution 
efficiently.  Currently,  when  a  satellite  is  added  to  the  navigation  solution,  the  MMAE 
is  completely  reinitialized  to  the  floating-point  filter,  and  new  candidate  ambiguity 
sets  are  generated.  The  removes  any  information  or  advantage  already  gained  by  the 
MMAE.  The  performance  of  the  algorithm  would  be  improved  by  the  using  the 
previous  MMAE  solution  to  reinitialize  all  of  the  elemental  filters  and  to  generate 
new  candidate  ambiguity  sets.  The  algorithm  might  also  be  changed  to  estimate  the 
ambiguity  of  the  new  measurement  accurately  and  incorporate  it  into  the  MMAE 
without  reinitializing  the  entire  MMAE. 

7.  Change  the  algorithm  to  a  moving  bank  MMAE  design.  Instead  of  fixing  the 
ambiguity  sets  hypothesized  by  the  elemental  filters,  allow  the  MMAE  to  move 
throughout  the  search  space  to  find  the  correct  elemental  filter. 

8.  Modify  the  current  algorithm  to  estimate/resolve  LI,  L2  and  widelane 
ambiguities  independently  and  simultaneously.  Widelane  ambiguities  are  a  linear 
combination  of  LI  and  L2  ambiguities,  and  the  additional  ambiguity  estimates 
could  be  used  to  validate  the  other  ambiguity  solutions,  and/or  increase  the 
accuracy  of  the  MMAE  final  solution. 
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Appendix  A.  Flight  Test  Profiles 


INTRODUCTION 

The  flight  tests  described  in  this  Appendix  were  flown  as  part  of  a  combined  Air 
Force  Institute  of  Technology  (AFTT)  and  USAF  Test  Pilot  School  (USAF/TPS)  research 
and  flight  test  program.  Flight  test  points  and  maneuvers  were  divided  into  primary  and 
secondary  flight  test  points.  The  specific  flight  test  points,  maneuvers,  and  techniques 
flown  during  this  thesis  are  described  in  the  following  sections.  Due  to  funding 
limitations  and  time  constraints  not  all  test  points  detailed  below  were  accomplished 
during  the  flight  test  phase  of  this  research.  Only  the  primary  flight  test  points  were 
accomplished. 

PRIMARY  FLIGHT  TEST  POINTS 

The  steady  flight  test  points  shown  in  Table  A-1  were  the  primary  data  points.  Figure 
A-1  gives  a  graphical  illustration  of  the  primary  flight  test  point  profiles.  These  flight  test 
points  started  over  the  reference  receiver  location  (Point  A),  flew  out  to  the  arc  distance 
(Point  B),  flew  the  desired  arc  (Point  C),  then  flew  a  straight  line  back  to  the  reference 
receiver  (Point  A).  The  arc  distance  was  defined  as  the  horizontal  distance  from  the  GPS 
reference  station  to  the  test  aircraft  (Point  A  to  B).  The  time  listed  for  each  test  point 
combination  of  arc  distance  and  altitude  was  the  estimated  total  time  required  to  fly  the 
entire  route  -  from  the  reference  receiver  to  the  arc,  the  arc,  and  then  back  to  the  reference 
receiver.  The  arc  time  was  the  estimated  length  of  time  data  will  be  collected,  at  the 
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required  arc  distance  and  altitude,  while  flying  each  arc.  The  arc  times  listed  in  Table 
A-1  were  based  on  a  0.6  Mach  number.  The  number  of  runs  listed  in  Table  A-1  was  the 


required  number  of  data  runs  at  each  test  point.  Primary  data  point  flight  test  tolerances 
are  listed  in  Table  A-2. 


Table  A-1.  Primary  Flight  Test  Data  Points 


Arc  Distance  | 

Pressure  Altitude 

10  nm 

20  nm 

30  nm 

60  nm 

100  nm 

10,000  ft 

■BSH 

X 

X 

20,000  ft 

■assM 

32  min 

X 

30,000  ft 

7  min 

15  min 

22  min 

32  min 

60  min 

Approx.  Arc  Time 

4  min 

8  min 

12  min 

12  min 

26  min 

Number  of  Runs 

3 

2 

2 

1 

1 

Table  A-2.  Primary  Flight  Test  Data  Point  Tolerances 


Parameter 

Tolerance 

Arc  Distance 

±  10%  of  Arc  Distance  or  5  nm  whichever  is  less 

Altitude 

±  1,000  ft 

Approx.  Arc  Time 

±  1  min  if  Arc  Distance  <30  nm,  otherwise  ±  2  min 

Figure  A-1.  Primary  Flight  Test  Point  Profile 


SECONDARY  FLIGHT  TEST  POINTS 

Djoiamic  flight  test  points,  or  secondary  data  points,  were  started  at  the  flight 
conditions  shown  in  Table  A-3.  The  letters  in  Table  A-3  at  each  test  point  combination 
of  arc  distance  and  altitude  refer  to  a  specific  flight  test  maneuver,  described  below. 
Secondary  data  point  tolerances  are  listed  in  Table  A-4. 

Lgve/  Turns 

A  constant  load  factor  level  turn  was  accomplished  to  investigate  angular  acceleration 
effects  on  position  and  velocity  accuracy.  These  turns  were  accomplished  at  load  factors 
of  1.4G,  2G  and  3G.  Each  maneuver  was  completed  at  an  airspeed  required  to  maintain 
the  constant  load  factor. 

Wind-Up  Turn 

A  constant  speed  Wind  Up  Turn  (WUT)  was  completed  to  investigate  angular 
acceleration  rate  effects  on  position  and  velocity  accuracy.  The  WUT  maneuver  varied 
load  factor  (G)  from  IG  to  a  maximum  of  6G  at  an  onset  rate  of  IG  per  4  seconds. 

Push  Over  /  Pull  Up 
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A  Push  Over  /  Pull  Up  was  completed  to  investigate  vertical  acceleration  effects  on 
position  and  velocity  accuracy.  The  Push  Over  /  Pull  Up  varied  load  factor  (G)  from  OG 
to  a  maximum  of  3G  at  an  onset  rate  of  IG  per  4  seconds. 


Table  A-3.  Secondary  Flight  Test  Data  Points 


Arc  Distance  | 

Pressure  Altitude 

10  nm 

20  nm 

30  nm 

10,000  ft 

A,B,C 

A,  B,C 

A,  B,C 

20,000  ft 

A,B,C 

A,B,C 

A,  B,C 

Table  A-4.  Secondary  Flight  Test  Data  Point  Tolerances 


Parameter 

Tolerance 

Arc  Distance 

Arc  distances  shown  (±  10%)  are  conditions  for  the 
start  of  the  maneuver.  The  Arc  Distance  may  vary  as 
necessary  to  complete  the  maneuver. 

Altitude 

±2,000  ft 

Airspeed  /  Mach 

±  10  KIAS/±  0.05  Mach 

Load  Factor 

±0.2G 

Descendins  Turn 

An  additional  secondary  flight  test  point  was  a  Descending  Turn.  This  maneuver 
investigated  altitude  effects  on  position  and  velocity  accuracy.  The  maneuver  was  flown 
at  an  arc  distance  of  5  (±1)  nm  from  the  reference  receiver,  and  a  VVI  of  4000  (±200) 
fpm.  The  maneuver  will  start  at  30,000  ft  PA  and  end  at  10,000  ft  PA. 

REFERENCE  RECEIVER  LOCATION 
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The  reference  receiver  station  for  all  flight  test  was  located  at  the  following  WGS-84 


coordinates 


Latitude 

(Degrees:Minutes:Seconds) 

Longitude 

(Degrees:Minutes:Seconds 

Altitude  (m) 

34:55:  19.88982 

-117:53:27.38584 

679.8870 
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